A method for generating dynamic compression shear coupled stress loading on living cells

Changes in the mechanical properties of single cells are related to the physiological state and fate of cells. The construction of cell constitutive equations is essential for understanding the material characteristics of single cells. With the help of atomic force microscopy, bio-image processing algorithms, and other technologies, research investigating the mechanical properties of cells during static/quasi-static processes has developed rapidly. A series of equivalent models, such as viscoelastic models, have been proposed to describe the constitutive behaviors of cells. The stress-strain relations under dynamic processes are essential to completing the constitutive equations of living cells. To explore the dynamic mechanical properties of cells, we propose a novel method to generate a controllable dynamical compression shear coupling stress on living cells. A CFD model was established to visualize this method and display the theories, as well as assess the scope of the application. As the requirements or limitations are met, researchers can adjust the details of this model according to their lab environment or experimental demands. This micro-flow channel-based method is a new tool for approaching the dynamic mechanical properties of cells.

Introduction

During the processes of cell development, differentiation, physiology, and disease, cells receive not only chemical signals but also mechanical signals from the extracellular matrix and surrounding environment (Hamed, 2020). Mechanical forces are experienced by cells and may be interpreted as a signal to induce phenotypic and functional responses or pathways, such as gene expression cascades, protein synthesis, proliferation, and movement; these responses can temporarily or even permanently change the cellular state (Desprat et al., 2005; Patterson et al., 2019). Moreover, the mechanobiological responses of biological cells had been extensively studied also, e.g., the responses of mesenchymal stem cells and chondrocytes to mechanical stimuli (Zhang et al., 2008; Zhang et al., 2012; Ganadhiepan et al., 2019). From a mechanical perspective, cells present a special material that is far more complex than ordinary materials, such as metal and glass. It is worth noting that the mechanical properties of cells remain unstable in most pathological processes, such as metastasis, asthma, sickle cell anemia, and apoptosis (Alberts et al., 2002; Desprat et al., 2005; Hao et al., 2020). Thus, understanding the mechanical behaviors of cells can provide a useful perspective on describing disease progression and revealing the fundamental mechanisms of the working of biomaterials (Bao and Suresh, 2003; Moeendarbary and Harris, 2014).

In exploring the mechanical behaviors of cells and establishing stress-strain relationships in them, it is a challenge to properly apply a controllable force on the tissue/monolayer/cell and capture its real-time strain at a single cell scale (i.e., ${10}^{-5}$ m) (Patterson, 2020). In this regard, various reasonable assumptions have been proposed, which present a research-scale perspective on mechanical methods, including the mechanical and biological methods (Moeendarbary and Harris, 2014; Hao et al., 2020).

In conjunction with atomic force microscopy (AFM), the bio-image processing algorithm (Dudani et al., 2013), micropipette aspiration (MA) (Evans and Yeung, 1989), and microfluidic platforms (Urbanska et al., 2020) are the most common and effective mechanical tools to apply compression/tensile or shear stress on cells. In addition, to improve accuracy and collect more information, some modified techniques and methods have been designed, such as magnetic twisting cytometry (MTC) (Trepat et al., 2004) and uniaxial stretching rheometer (USR) (Desprat et al., 2005). Through these static or quasi-static mechanical experiments, it is believed that the mechanical behavior of cells is likely to resemble that of viscoelastic material (Katti et al., 2019). However, experiments exploring the stress-strain relationship of isolated living cells did not reach dynamic conditions or higher strain rates ($\dot{ϵ}>{10}^{-1}$ s⁻¹). Commonly, the dynamic loading process of materials, including living cells, differs significantly from the static or quasi-static situations. For instance, a quasi-static deformation situation comprises a sequence of equilibrium states, where the well-known equations describing the mechanics of materials work (i.e., $\sum F=0$; $\sum M=0$ ). On the contrary, the dynamic loading process can be treated as a stress wave traveling through the body at an acoustic speed (Meyers, 1994). When an external deformation is imparted at a very high rate, it induces stress on one portion of the body, while other portions remain unaffected. As cells can sense mechanical behaviors, they react rapidly to adapt to them (Pelham and Wang, 1997). Cells subjected to dynamic loading processes exhibit distinct mechanical behaviors rather than viscoelastic material. Moreover, the stress-strain relations under dynamic processes are an important part of the single-cells constitutive equations. Therefore, developing methods to apply dynamic stress on cells will significantly promote the understanding of the mechanical behaviors of cells under dynamic conditions (Bao and Suresh, 2003).

In the present work, we propose a novel method to apply combined dynamic compression-shear loading on isolated living cells under normal conditions; in addition, we extend the range of stiffness tensor of single biological cell piecewise function to higher strain rates ($\dot{ϵ}>{10}^0$ s⁻¹). The basis of this new method is the theory of transient flow, or in detail, the theory of the weak shock wave (where “weak” signifies that the thermal energy generated by impact compression is smaller compared to the total internal energy of the fluid (Smith, 1973)) propagation in a viscous fluid, which would suddenly induce both the compression and shear stresses in the boundary layer. In addition, the water hammer theory asserts that we can repeatedly apply disturbance with the amplitude of the weak shock wave, which is precisely controllable by changing the speed of the projectile. Briefly, this method includes two parts: the stress loading part accelerates a projectile to impact a fluid-fulfilled microchannel seeded with living cells on the bottom, while the strain acquisition part is equipped with a high-speed camera to assist with an image processing algorithm. Once the assumptions and requirements are fulfilled, the details of the design remain readily changeable. The proposed methods can be used to explore the stress-strain relations under dynamic processes and clarify the constitutive behavior of single cells to dynamic loadings.

Methods

To illustrate the dynamic and coupled compression-shear loading method, a simple schematic diagram of the model is presented in Figures 1A,B. The model exhibits a gas gun, projectile, and a microfluidic chip with a rectangular conduit channel, where cells can be seeded on the bottom wall. In this mode, a projectile was accelerated to an initial velocity ( $u_p$ ) to impact the buffer; a pressure surge was then generated, which traveled through the fluid matter. As this stress wave propagated at the acoustic velocity, the cells cultured on the bottom wall experienced the compression stress by the stress wave directly, as well as shear stress due to the viscosity of the fluid (Figure 1C).

FIGURE 1

FIGURE 1. A schematic diagram showing the basics of the proposed model. (A) The schematic diagram of the geometry of the model; (B) The sketch of the rectangular channel; (C) The stress analysis of cells during wave propagation.

Navier-Stokes equations

To allow the pressure waves to propagate, the fluid used in this experiment was compressible with constant viscosity. Therefore, the equations of continuity and momentum (Navier-Stokes equations) used to describe the motion of the fluid are written as:

$\begin{array}{l}\frac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho u\right)=0\end{array}(1)$

$\begin{array}{l}\frac{\partial \left(\rho u\right)}{\partial t}+\nabla \cdot \left(\rho u\cdot u\right)=-\nabla P+\mu {\nabla }^{2}u+\frac{\mu }{3}\nabla \left(\nabla \cdot u\right)\end{array}(2)$

Where $\rho$ is the fluid density; $\mu$ is the dynamic viscosity; $u$ is the fluid velocity; $″\nabla ″$ is the gradient operator; $″\nabla \cdot ″$ is the divergence operator; $″{\nabla }^2″$ is the Laplace operator which means $″\nabla \cdot \nabla ″$ ; $t$ is time, and $P$ is the total pressure in the fluid. For isotropic and homogeneous newton fluid, the stress tensor can be expressed as:

$\begin{array}{l}{\sigma }_{ii}=-P+2\mu \frac{\partial u_i}{\partial x_i}-\frac{2}{3}\mu \nabla \cdot u\end{array}(3)$

$\begin{array}{l}{\tau }_{ij}=\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)\end{array}(4)$

Where $u_i$ is the fluid velocity, and $x_i$ is the spatial coordinate ($i=1,2,3$). ${\sigma }_{ii}$ is the pressure components of stress tensor, and ${\tau }_{ij}$ is the shear components of stress tensor. In the channel flow with non-slip wall condition, the wall shear stress $\left({\tau }_i^w\right)$ could be calculated using:

$\begin{array}{l}{\tau }_i^w=\mu \left(\frac{\partial u_i}{\partial n_w}\right)\end{array}(5)$

Where $n_w$ is the distance normal to the wall.

Wave propagation celerity in a rectangular conduit channel

A general solution for the celerity of pressure wave propagating in the fluid has been derived by Hutarew and can be applied for a conduit of any cross-section (Hutarew, 1973). The form is:

$\begin{array}{l}{c}_f=1/\left(\rho \left(\frac{1}{K}+\frac{1}{\mathcal{A}}\frac{\delta \mathcal{A}}{\delta p}\right)\right)\end{array}(6)$

Where K is the fluid bulk modulus, $\mathcal{A}$ is the cross-section area of the channel, and $c_f$ is the wave propagation celerity. $\delta$ is a small change sign. Eq. 6 shows the effect of fluid-structure interaction on the propagation speed of the pressure surge. This means that despite the fluid in the Computational Fluid Dynamics (CFD) model being considered compressible, the deformation of channel walls should be taken into account. However, in this study, a rigid-walled boundary condition was preferred, which could make the measurement and calculation of the cell strain field much more concise. Hence, a corrected fluid bulk modulus $\left(K_{\mathit{mod}}\right)$ was proposed to incorporate structural behavior in he CFD model. For the propagation of pressure surge in liquids traveling through thick-walled pipes and ducts of rectangular cross-section, the theoretical wave celerity equation was given by Thorley (Thorley and Guymer, 1976):

$\begin{array}{l}{c}_f=1/\left(\rho \left(\frac{1}{K}+\frac{\Phi \left(a,b\right)}{abE{e}^2}\right)\right)\end{array}(7)$

$\Phi \left(a,b\right)=\frac{a^3+b^3}{2\left(a+b\right)}\left(\frac{a^3}{6}+\frac{a^{2}b}{2}-\frac{b^3}{3}\right)-\frac{a^5}{20}\begin{array}{l}-\frac{a^2{a}^3}{4}+\frac{b^5}{5}+\frac{E{e}^2}{4G}\left(a^3+b^3\right)+\frac{ab{e}^2}{2}\left(a+b\right)\end{array}(8)$

Where a and b are the lengths of the long and short sides of the rectangular cross-section; e is the wall thickness; $E$ and G are the elastic and shear moduli of the wall material, respectively. $\Phi \left(a,b\right)$ was solved using Eq. 8. According to the equation defining wave propagation velocity, the corrected fluid bulk modulus $K_{\mathit{mod}}$ equation is:

$\begin{array}{l}{c}_f={\left(\frac{K_{\mathit{mod}}}{{\rho }_0}\right)}^{\frac{1}{2}}\end{array}(9)$

Nevertheless, the properties of solids involved played an important role in wave propagation as per the wave celerity equation. An evaluation of the fluid-solid coupling effect was needed to examine the limitations when the motion of the channel was rigid. A dimensional parameter “$\beta$”, named “fluid loading”, was proposed by Pinnington to evaluate the fluid-solid coupling; it corresponded to the Korteweg wave speed equation (Korteweg, 1878; Pinnington, 1997; Shepherd and Inaba, 2009).

$\begin{array}{l}\beta =\left(c_0/c_f\right)\end{array}(10)$

$\begin{array}{l}\beta =\left(\frac{c_f^2}{c_s^2}\right)\left(\frac{{\rho }_f}{{\rho }_s}\right)\left(\frac{2R}{e}\right)\end{array}(11)$

Where $c_0$ is the acoustic speed in fluid; $R$ is the mean radius of the channel; the subscript, “$s$”, denotes a structure, and the subscript, “$f$”, means fluid. Moreover, according to the limiting cases of fluid-structure interaction (FSI) discussed by Shepherd and Inaba, the case, where $\ll 1$ , indicated that the channel could be regarded as rigid (Shepherd and Inaba, 2009). Substituting the wave celerity and acoustic speed in Eq. 10 to get the parameter $\beta$ helped us design the details of the channel.

Pressure wave generation by projectile impact

To describe the fluid suitable here with equations, the Tait and Tammann equations of state, which apply to a wide range of liquids, were used (Dymond and Malhotra, 1988; Hoang et al., 2015). The Tait equation can be written as a function of pressure $\left(p\right)$ and density $\left(\rho \right)$

$\begin{array}{l}p=B{\left(\frac{\rho }{{\rho }_0}\right)}^{\gamma }-B\end{array}(12)$

Where $B$ is a weak function of entropy, which is usually treated as a pressure constant of $3.35\times {10}^8$ Pa. $\gamma$ is the adiabatic exponent that equals $7.15$ (Nagayama et al., 2002). The Tait equation was rewritten as a partial differential equation form.

$\ln \left(p+B\right)=\ln \left(B\right)+\gamma \ln \left(\frac{\rho }{{\rho }_0}\right)$

$\begin{array}{l}\frac{\partial p}{\partial \rho }=\frac{\gamma \left(p+B\right)}{\rho }\end{array}(13)$

As the whole process was assumed to be isentropic, the adiabatic exponent $\gamma$ was:

$\begin{array}{l}\gamma =-\frac{\gamma }{p}{\frac{dp}{d\gamma }|}_s\end{array}(14)$

Where γ is the volume Upon substituting volume γ by density ρ, we obtained:

$\begin{array}{l}\frac{p}{p_0}={\left(\frac{\rho }{{\rho }_0}\right)}^{-\gamma }\end{array}(15)$

The definition of acoustic velocity is:

$\begin{array}{l}{c}^2=\left(\frac{\partial p}{\partial \rho }\right)\end{array}(16)$

Upon integration of Eq. 16, a function of pressure and wave speed describing the fluid density was obtained, which could be used to correct the density in CFD, as the fluid was considered compressible.

$\begin{array}{l}\rho ={\rho }_0+c^{-2}p\end{array}(17)$

The acoustic velocity is derived from Eq. 15.

$\begin{array}{l}{\left(\frac{\partial p}{\partial \rho }\right)}_s=\gamma \left(\frac{p_0}{{\rho }_0}\right){\left(\frac{\rho }{{\rho }_0}\right)}^{\gamma -1}=c^2\end{array}(18)$

As the compressibility of the fluid was limited, an approximation of the acoustic velocity under normal conditions could be written as:

$\begin{array}{l}\underset{\rho \to {\rho }_0}{\lim }{\left(\frac{\partial p}{\partial \rho }\right)}_s=\frac{\gamma p_0}{{\rho }_0}=c_0^2\end{array}(19)$

The local acoustic velocity could be written as:

$\begin{array}{l}c={\left(\frac{\partial p}{\partial \rho }\right)}^{\frac{1}{2}}=c_0{\left(\frac{p}{p_0}\right)}^{\frac{\gamma -1}{2\gamma }}=c_0{\left(\frac{\rho }{{\rho }_0}\right)}^{\frac{\gamma -1}{2\gamma }}\end{array}(20)$

Figure 2A displays the shock wave jumps conditions with a coordinate fixed at the shock wave front. The equations of conservation (mass; momentum; energy) could be written as:

$\begin{array}{l}{\rho }_0\left(u_0-c_s\right)={\rho }_1\left(u_1-c_s\right)\end{array}(21)$

$\begin{array}{l}{\rho }_0\left(u_0-c_s\right)u_{0}dt={\rho }_1\left(u_1-c_s\right)u_{1}dt\end{array}(22)$

$\begin{array}{l}\frac{1}{2}{u}_0^2+\frac{\gamma }{\gamma -1}\frac{p_0}{{\rho }_0}=\frac{1}{2}{u}_1^2+\frac{\gamma }{\gamma -1}\frac{p_1}{2}\end{array}(23)$

FIGURE 2

FIGURE 2. The isolation of boundary conditions and analysis of the model. (A) The weak shock wave jumps conditions under a coordinate fixed with the shock wave front; (B) The weak shock wave generated by projectile impact analysis; (C) A schematic diagram illustrating the micro-channel boundary conditions in CFD modeling.

These equations are usually known as the Rankine-Hugoniot equations (Smith, 1973). Here, we neglected the thickness of the shock front, and thus, the particle velocity, $u$ and shock velocity, $c_s$ behind the shock front could be solved. For very weak shock waves, the jump process could be treated as an isentropic process, yielding the following equation:

$\begin{array}{l}{c}_s={\left(\left(\frac{\rho }{{\rho }_0}\right)\left(\frac{p-p_0}{\rho -{\rho }_0}\right)\right)}^{\frac{1}{2}}\end{array}(24)$

$\begin{array}{l}u=\frac{\left(\rho -{\rho }_0\right)}{\rho }{c}_s\end{array}(25)$

$\begin{array}{l}\frac{{\rho }_1}{{\rho }_0}=\frac{\left(\gamma +1\right)p_1+\left(\gamma -1\right)p_0}{\left(\gamma +1\right)p_0+\left(\gamma -1\right)p_1}=\frac{u_0}{u_1}\end{array}(26)$

The particle velocity could then be rewritten with Eq. 19.

$\begin{array}{l}u=\frac{p-p_0}{{\rho }_0{c}_0}=\frac{c_0}{\gamma p_0}\left(p-p_0\right)\end{array}(27)$

or

$\begin{array}{l}\frac{p}{p_0}=1+\frac{\gamma u}{c_0}\end{array}(28)$

Substituting by Eq. 20, gave:

$\begin{array}{l}{c}_s=c_0{\left(1+\frac{\gamma u}{c_0}\right)}^{\frac{\gamma -1}{2\gamma }}\end{array}(29)$

The Taylor expanded form of Eq. 29 thus obtained is:

$c_s=c_0\left(1+\frac{\gamma -1}{2c_0}u+\mathcal{O}\left(u^2\right)+\dots \right)\cong c_0\left(1+\frac{\gamma -1}{2c_0}u\right)\begin{array}{l}\cong c_o+\left(\frac{\gamma -1}{2}\right)u\end{array}(30)$

As the relationship between acoustic velocity and particle velocity was obtained, the differential form was:

$\begin{array}{l}d{u}_p=du=\frac{2}{\gamma -1}d{c}_S\end{array}(31)$

As no cavitation or fluid column separation occurred and no cross-section changes were recorded along the channel, the Joukowsky equation was a perfect approximation to predict the maximum pressure in a water-hammer impact situation (Joukowsky, 1900; Walters and Leishear, 2018). This impact occurs at $\left(x,t\right)=\left(\mathrm{0,0}\right)$, and an equation, $p\left(t\right)$, describes the dynamic pressure of the fluid at the interface of the fluid and projectile. This initial impact would create a weak shock wave with an amplitude that can be determined by the pressure-velocity matching method (Meyers, 1994; Deshpande et al., 2006; Shepherd and Inaba, 2009)

$\begin{array}{l}\Delta p=p\left(0\right)=\frac{{\left(\rho c\right)}_f{\left(\rho c\right)}_p}{{\left(\rho c\right)}_f+{\left(\rho c\right)}_p}{u}_p\end{array}(32)$

Where ${\left(\rho c\right)}_f$ is the fluid impedance; ${\left(\rho c\right)}_p$ is the projectile impedance. The impedance of a projectile is much higher than that of the fluid in most cases, and Eq. 12 can reduce it to an approximated expression:

$\begin{array}{l}p\left(0\right)\approx {\left(\rho c\right)}_f{u}_p\end{array}(33)$

After the conditions of the fluid were solved, the impact condition obtained is shown in Figure 2B. The shock wave generation due to the projectile impact was a discrete process. In this study, since the buffer was assumed to be as thin as possible, it could be neglected given that it had the same impendence as the projectile. This model could be simplified to a projectile impacting the fluid column directly. Analysis of the projectile motion was performed according to Newton’s second law.

$\begin{array}{l}{p}_i{\mathcal{A}}_p=m_p\frac{d{u}_p}{dt}\end{array}(34)$

Where $p_i$ is the pressure by every discrete impact; $u_p$ is the interface velocity of projectile and fluid; ${\mathcal{A}}_p$ is the cross area of the projectile, and $m_p$ is the mass of the projectile. Substituting Eq. 34 with Eq. 31, the result was:

$\begin{array}{l}p=\frac{2}{\gamma -1}\frac{m_p}{{\mathcal{A}}_p}\frac{d{c}_s}{dt}\end{array}(35)$

Meanwhile, Eq. 12 was rewritten in terms of local acoustic velocity $c_s$ to obtain the following:

$\begin{array}{l}p=B{\left(\frac{{\rho }_0{c}^2}{nB}\right)}^{\frac{\gamma }{\gamma -1}}-B\end{array}(36)$

The partial differential form of Eq. 33 is:

$\begin{array}{l}\frac{\partial p}{\partial c}=\frac{2\gamma }{\gamma -1}\frac{P+B}{c}\end{array}(37)$

In our model, the pressure was too small as compared with the pressure constant, which is $3.35\times {10}^8$ Pa (Nagayama et al., 2002). A linear approximation treatment (Lennon, 1994) was used with normal condition parameters $\left(p_0,c_0\right)$ to calculate Eq. 34, which could be rewritten as:

$\begin{array}{l}{\left(\frac{\partial p}{\partial c}\right)}_s=\frac{2\gamma }{\gamma -1}\frac{P_0+B}{c_0}\end{array}(38)$

The normal conditions here were $p_0=101325Pa$ , ${\rho }_0=999.8$ kg/m³, and $c_0=1439$ m/s. Eq. 33 was then rewritten as follows:

$\begin{array}{l}p=\frac{2}{\gamma -1}\frac{m_p}{{\mathcal{A}}_p}\frac{d{c}_s}{dt}=\frac{2}{\gamma -1}\frac{m_p}{{\mathcal{A}}_p}{\left(\frac{\partial c}{\partial p}\right)}_s\frac{dp}{dt}\end{array}(39)$

Upon integration of Eq. 36, a simplified function of pressure was gained.

$\begin{array}{l}p\left(t\right)=p\left(0\right)\exp \left(-\lambda t\right)\end{array}(40)$

Where the time constant was as follows:

$\begin{array}{l}\lambda =\frac{\gamma -1}{2}\frac{{\mathcal{A}}_p}{m_p}{\left(\frac{\partial p}{\partial c}\right)}_s\end{array}(41)$

Boundary layer induced by a very weak shock wave

Once a weak shock wave enters a static fluid boundary through the wall, a boundary layer begins to appear at the interaction point of the weak shock wave and the wall (Ackroyd, 1967; Davies and Bernstein, 1969). A study of the laminar compressible boundary layer induced by this weak shock wave was solved by H. Mirels’ in 1955, who provided theoretical evidence to prove that the weak shock waves generated by projectile impact indeed induce a boundary layer (Mirels, 1955).

Considering a possible turbulent boundary layer, R. E Melnik and B. Grossman developed an asymptotic theory within the limit of weak shock waves (Melnik and Grossman, 1974). Their three-layer description of the boundary layer was a natural extension of the asymptotic theories of Mellor (Mellor, 1972), Yajnik (Yajnik, 1970), Bush and Fendell (Bush and Fendell, 1972; Bush and Fendell, 1973) for incompressible boundary layers, as well as the theory by Afzal (Afzal, 1973) for compressible non-interacting turbulent boundary layers. FLUENT is a reliable and acceptable CFD software employed for relevant numerical simulations. To simulate our model, the k-ε model was utilized to analyze turbulent flow (Versteeg and Malalasekera, 2007; Nikpour et al., 2014).

Geometry definition

In the present work, the geometry definition consisted of a volume occupied by the flow with a specified shape of the physical boundaries. Here are several important criteria to limit the size of the channel, although, in the beginning, we planned to put the ‘lab’ into just one small chip to make this as convenient as possible. A fundamental restriction applied was the continuous flow in the fluid domain to ensure that the theories being considered were effective. A dimensionless parameter, the Knudsen number, was used to describe this problem (Guo et al., 2013); this number is defined as

$\begin{array}{l}{K}_n=\mathrm{ƛ}/{\mathcal{\ell }}_{char}\end{array}(42)$

Where $\mathrm{ƛ}$ is the mean free path of the particle, and ${\mathcal{\ell }}_{char}$ is the characteristic length scale. The fluid domain was described as a continuum and solved by Navier-Stokes (N-S) equations with no-slip boundary conditions, therefore, the $K_n$ was generally considered to be less than 0.001 (Ben-Dor et al., 2000; Kandlikar et al., 2005; Rapp, 2016).

According to the requirements above, we established a simple CFD model to describe our method. The Cartesian coordinates system was created as shown in Figure 1B, and a $3\times 3\times 60$ mm rectangular conduit channel was drawn in ANSYS. The details of size were kept adjustable to adapt to various laboratory environments for other researchers. An example size is proposed here to describe dynamic stress loading methods.

Mesh generation and solver settings

Transient simulation is strongly dependent on the quality of the mesh. For most water-hammer models to accurately capture the near-wall velocity, the mesh density near the wall should be concentrated. During the cross-section meshing process, the boundary layer neighboring the wall was divided into 20 layers with the first layer having a thickness of $1\times {10}^{-6}$ m (1.1 growth rate); this design was based on the mesh independence analysis of 3D pressurized pipe flow with CFD modeling previously described by Martins et al. (Martins et al., 2018). Given the axial direction, a sweep method was applied. The Courant-Friedrichs-Lewy (CFL) criterion was used to maintain stability during the movement of the acoustic wave across the discrete elements in CFD, which determined the element size in the axial direction. The non-dimensional Courant number was calculated using the following equation:

$\begin{array}{l}{\mathit{C}}_{\mathrm{c}\mathrm{o}\mathrm{u}\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{t}}=c_f\frac{\Delta t}{\Delta z}\end{array}(43)$

Where $\Delta z$ is the element size in the axial direction, and $\Delta t$ is the calculation time step. There were two main considerations for the determination of $C_{courant}$: one was to capture the wave velocity in the fluid; and the other was to calculate the stress on an arbitrary cell at the channel wall. The $C_{courant}$ was ideally expected to be ≤ 1. Moreover, to maintain the meshing quality, i.e., to keep the aspect ratio of all mesh elements < ${10}^3$, the axial element size was $\Delta z=0.1$ mm with a time step of $t=0.7\times {10}^{-7}$ s. The total time was assumed to be $4.2\times {10}^{-5}$ s, i.e., 600 steps in total, to ensure that the entire compression wave could travel completely from inlet to outlet. The total mesh count was 1.98 million elements (Martins et al., 2014; Mandair, 2020).

In the CFD solver setup, boundary settings comprise the physical description of the fluid flow as shown in Figure 2C. The pressure applied on the inlet boundary has been described as a function of time, $\left(t\right)$ in Eq. 40, which simulated the impact process. The outlet was operated under normal pressure ($p_0=10135$ Pa), which could be achieved easily by collecting an extra tank filled with fluid, and the wall was set at a no-slip condition. In addition, there were two settings (heights) for the roughness of the bottom wall: $1\times {10}^{-5}$ (average height of cells) and $0$ m, which helped estimate the effect of cells on the fluid flow.

ANSYS Fluent^® (2019R3) was utilized to obtain all the simulation results presented in this paper. In FLUENT, the numerical tech is a finite volume method (FVM). The whole process is transient. As the fluid is considered to be viscous, compressible, isothermal (no heat transfer), isotropic, and single-phase (no cavitation), the semi-implicit method for pressure-linked equations (SIMPLE) can be used as a flow solver. Convergence was defined to be $1\times {10}^{-6}$ due to the flow characteristics (Martins et al., 2016).

Results

Sample numerical results

A mesh independence analyses were performed using different element size (the original element size was proposed in Section 2.6) by the simulation of pressure and shear stress on bottom wall. The total mesh count for testing the independence were about $0.855\times {10}^6$ (blue line), $1.98\times {10}^6$ (green line), $10.134\times {10}^6$ (red line) elements respectively. The test results were shown in Figures 3A,B which assure a mesh independence of the simulation results in this paper.

FIGURE 3

FIGURE 3. Mesh independence analysis and sample comparison results. The total mesh count for testing the independence were about $0.855\times {10}^6$ (blue line), $1.98\times {10}^6$ (green line), $10.134\times {10}^6$ (red line) elements respectively. (A) The variations of total pressure on bottom wall at the monitoring point ( $x=0.25X$); (B) The variations of axial shear stress on bottom wall at the monitoring point ( $x=0.25X$); (C) A pressure waveform comparison of experimental results of example and numerical results (black line: experimental results; blue line: numerical results).

We choose a classical water hammer experiment (Inaba and Shepherd, 2010) to assess the modeling method in our work. In this example (shot 62), researchers accelerated a steel projectile ($=0.67$ kg, $v_{impact}=18.5$ m/s) to impact a specimen tube ($D_{innerdiameter}=38.1$ mm, $h_{thickwall}=12.74$ mm), and the experimental pressure was recorded by strain gauges and a pressure transducer.

In Figure 3C, we compared the experimental results (black line; shot 62, gauge g1) with the numerical results (blue line) calculated with the modeling method in Section 2. The maximum pressure of $25.03$ MPa was computed which showed a good match with the peak pressure of $27.20$ MPa. Furthermore, the whole trend of the pressure waveform showed a good agreement with the experimental results. For instance, due to the outlet of the specimen tube being closed, a reflected wave could be observed.

Visualization of pressure wave propagation

To visualize and analyze the wave propagation process in detail, as well as the profile of the coupled compression and shear stresses in our method, a representative case with a water-filled channel and a projectile made of polymethyl methacrylate (PMMA) was designed; a practical size was determined as described in Section 2.5. The correlation between the projectile parameters (velocity and length) and the input pressure was assessed, and the results are shown in Figure 4. The maximum value and the duration of the pressure input were adjustable by changing the length and initial velocity of the projectile. This detailed relationship was derived in Section 2.3 (Joukowsky, 1900; Versteeg and Malalasekera, 2007).

FIGURE 4

FIGURE 4. Evaluation of the relationship between various lengths/velocities of the projectile and the input pressure profile. (A) Three different lengths of projectile with a 1 m/s velocity (green line: 4 cm; blue line: 2 cm; red line: 1 cm); (B) Three different lengths of projectile with a 1 cm length (green line: 4 m/s; blue line: 2 m/s; red line: 1 m/s).

As the inlet boundary condition, maximum pressure of $3$ MPa was generated by a $0.01m$ projectile, which was then used to initiate the wave propagation process. Along the direction of the wave propagation ($z$ -direction), a series of axial middle cross-sections were selected to visualize the wave traveling process in terms of pressure and axial flow velocity distribution at different time steps (Figure 6). Here, we set $={\mathcal{l}}_{tube}/c_f$, and the displayed time steps were $t_1=0.25T,t_2=0.50T,and{t}_3=0.75T$. After the impact of the projectile, the pressure jumped to ∼ $3.2$ MPa rapidly, accompanied by a relatively low level of flow velocity (less than $2.3$ m/s). As shown in Figure 5, peak pressures at the different time steps did not show an apparent dissipation, while the maximum axial flow velocity slow down slightly. Details on the dissipation values are discussed in Section 3.3. The phase differences calculated by Eq. 7 showed that the wave traveled at an acoustic speed.

FIGURE 5

FIGURE 5. A CFD visualization of weak shock wave propagation in the middle axial cross-section at different time steps ( $t_1=0.25T,t_2=0.50T,t_3=0.75T$ ). (A) The total pressure field; (B) The axial velocity field.

Stress distribution analysis of wall/cell cultured area

The bottom wall was defined as a smooth wall area for cell culture. The results at $t=0.75T$ were selected to display the stress distribution in space. Figure 6A shows the compression pressure applied on the bottom wall, which appears to be similar to the compression pressure distributed in the axial plane in Figure 5A. The pressure distributed was almost identical to that in the transverse direction ( -direction) (Figure 6B). The axial shear stress (${\tau }_z$) is shown in Figure 6C, which exhibited no obvious difference in the transverse direction. The maximum axial shear stress at $t=0.75T$ was nearly $1.9$ kPa, which was much higher than the maximum transverse shear stress ($\left|{\tau }_{x,\mathit{max}}\right|<1$ Pa). From these observations, it was inferred that the transverse shear stress (${\tau }_x$) could be neglected, even though it showed an apparent concentration of stress at the edges (Figure 6B).

FIGURE 6

FIGURE 6. The CFD visualization of stress distribution on the bottom wall where cells were cultured at $t=0.75T$ time step. (A) The total compression stress distribution; (B) The transverse shear stress distribution; (C) The axial shear stress distribution.

Effect of cells on stress spatial distribution

To evaluate the influence of cultured cells on the spatial distribution of stress on the bottom wall, the relationship between the two distributions and the density and height of the cells could be written as:

$\begin{array}{l}\sigma =f\left({\delta }_{cell},{\rho }_{cell}\right)\end{array}(44)$

Where $\sigma$ is the compression and shear stress on the bottom wall; ${\delta }_{cell}$ is the average height of the cells, and ${\rho }_{cell}$ is the density of the cells. Normally, the height of the cells, ${\delta }_{cell},$ was consistent at ∼ $1\times {10}^{-5}$ m and the interval of variation, ${\rho }_{cell},$ was $[0,\infty )$. The two boundary conditions were defined as follows: the lower limit situation, i.e., ${\rho }_{cell}=0$ (no cell cultured), represented a smooth wall; the upper limit situation, i.e., ${\rho }_{cell}=\infty$, denoted that countless cell had been cultured on the bottom wall and imparted a roughness of $1\times {10}^{-5}$ m. By comparing the stress in these two situations, the effect of cells on stress distribution could be assessed (Figure 7). Along the axial and transverse directions, several monitoring points were set to calculate the compression and shear stress: $\sigma$ (3 points along the axial direction at $x=0.5X$: $z_1=0.25\mathcal{l}$, $z_2=0.5\mathcal{l}$, $z_3=0.75\mathcal{l}$; 3 points along the $x$ direction at $z=0.5L$: $x_1=0.25X$, $x_2=0.5X$, $x_3=0.75X$).

FIGURE 7

FIGURE 7. The comparison of rough and smooth model historical variations in stress were recorded by monitoring the points calculated by FLUENT (blue line: rough model; red line: smooth model). (A) The transverse shear stress recorded at three monitoring points ($x_1=0.25X$, $x_2=0.5X,$ $x_3=0.75X$) at the line $z=0.5L$; (B) The axial shear stress recorded at three monitoring points ($z_1=0.25\mathcal{l}$, $z_2=0.5\mathcal{l}$, $z_3=0.75\mathcal{l}$) at the line $x=0.5X$: (C) Compression stress recorded at three monitoring points ($z_1=0.25\mathcal{l}$, $z_2=0.5\mathcal{l}$, $z_3=0.75\mathcal{l}$) at the line $x=0.5X$.

In addition, the temporal variations in ${\tau }_x$ in the transverse direction were also recorded. As the maximum of ${\tau }_x$ was quite small (did not exceed ${10}^1$ Pa) within the interval of $\left(0.25X,0.75X\right)$ , ${\tau }_x$ could be neglected (Figure 7A). As shown in Figures 7B,C, the compression and axial shear stress variations along the axial direction were recorded at the monitoring points, and the same changing trend was observed. Exact statistics on the maximum values are shown in Table 1. Regardless of the maximum value of pressure or axial shear, there was almost no difference between these two conditions in an interval of $\left(0.25\mathcal{l},0.75\mathcal{l}\right)$. The ${\rho }_{\mathit{exp}}$ was in the interval of $0<{\rho }_{\mathit{min}}\le {\rho }_{\mathit{exp}}\le {\rho }_{\mathit{max}}\ll \infty$. Therefore, the smooth boundary condition would be much closer to the actual experimental observations.

TABLE 1

TABLE 1. Maximum values of the monitoring points were calculated by FUENT.

Table 1 shows the attenuation of pressure and axial shear stress along the axial direction. In this model, the average reduction was within 1% in every $0.25\mathcal{l}$ ($15$ mm) for pressure and ∼ $6\%$ in every $0.25\mathcal{l}$ for shear stress.

In addition, the temporal variations in the transverse direction are shown in Figure 8, where the monitoring points were set in the middle cross-section ($z=0.5\mathcal{l}$). The attenuation of either the pressure or the axial shear stress in the transverse direction was less than $1\%$ . These results imply the great repeatability of this model.

FIGURE 8

FIGURE 8. The smooth model historical variations in stress recorded by monitoring points calculated by FLUENT (blue line: rough model; red line: smooth model). (A) The axial shear stress recorded at three monitoring points ( $x_1=0.25X$, $x_2=0.5X,$ $x_3=0.75X$) at the line $z=0.5L$; (B) The compression stress recordings at three monitoring points ($x_1=0.25X$, $x_2=0.5X,$ $x_3=0.75X$) at the line $z=0.5L$.

Discussion

In the past, many studies have focused on revealing the mechanical properties of cells with the assistance of interdisciplinary techniques, such as solid mechanics, image processing, and fluid mechanics. Such works have implied that as compared to static or quasi-static loadings, the mechanical properties of cells are likely to show a distinct behavior difference under dynamic loading ($\dot{ϵ}>{10}^0{s}^{-1}$). Considering the length scale of cells, it is difficult to clamp a single cell using traditional dynamic loading methods, such as Split Hopkinson’s Pressure Bar (SHPB). Hence, the challenge is how to apply controllable dynamic loading on cells properly and the real-time measurement of cell response during this process.

In this paper, we have proposed a new method for the application of dynamic loading that couple’s compression and shear stress synchronously on isolated cells under normal culture conditions. This method describes a case, where an impact loading from a projectile was employed to generate a pressure wave, together with corresponding shear stress due to the fluid viscosity. At the same time, the strain variations in the cell could be captured by a high-speed camera. Eventually, data on the essential stress and strain on the cells were collected to explore the mechanical properties of cells. We established a representative model (rectangular channel filled with water) as shown in Figures 1A,B. With this model, we explained two main questions as elaborated on in Section 2 and Section 3 on the working of this method and its known limitations:

1) How can the weak shock wave generated by projectile impact be calculated? Does the solid structure (channel) affect the flow and how can it reduce the fluid-structure interaction (FSI)?

2) Does the shear stress keep the same phase as the compression stress? Does the existence of cells have a strong influence on stress distribution? How is this related to the positions along the axial and transverse directions?

The basic purpose was to explore the biomechanical mechanism of single cell response to dynamic mechanical loading, where the length scale was focused on ${10}^{-6}$ m and the cells were treated as a homogenous element. This length scale had to be taken into consideration due to its influence on the fluid domain meshing and CFD (actually, more factors had been considered during the meshing process). The characteristic length of the channel we designed was $3\times {10}^{-3}$ m, which was much larger than the size of the cells. Generally, the pressure and wall shear stress were directly assumed when the stresses were applied to the cells, while the approximate treatment relied on the hypothesis that the flow was hydraulically smooth. The estimated average height of the cells was ∼ $1\times {10}^{-5}$ m, therefore, one of the extremely idealized hypotheses was that the bottom wall was occupied by countless cells that were equivalent to a $1\times {10}^{-5}$ m roughness on the bottom wall. The total pressure and axial shear stress showed that there was only a slight difference between the smooth and rough conditions (within $2\%$) (Figure 7). In the actual experimental process, the cell density only sparsely facilitated the capture of the strain field by a high-speed camera. Therefore, the effect of cells on the fluid flow could be neglected.

The Doppler effect is commonly used to detect the flow velocity of a flowing fluid, but experimental techniques may render it less reliable on small velocities and the changing profile near the walls (Riasi et al., 2009). Several sophisticated numerical models have been established to more accurately explain the details of the transient flow (Ghidaoui and Kolyshkin, 2001; Martins et al., 2018). It has been proven that CFD performs very well in modeling the pressure wave traveling processes (Mandair, 2020). Therefore, we established the CFD model to help us evaluate the proposed method, and we will consider these evaluation results as a very important reference for future experimental work.

The boundary conditions should be comprehensively considered in CFD modeling as we have described in detail in Section 2 of this article. The criteria of channel size design were not too strict, which allowed us to adjust the details freely according to the actual lab environment, provided that the size details obeyed the above-mentioned limitations of the Knudson number to continuously maintain the fluid. As for the projectile, its material should be kept the same as that of the buffer. Figure 4 shows the means of controlling the amplitude and decay time constant of the initial pressure by the projectile.

The water hammer is a well-known problem. The weak shock wave generated in this problem would induce a sudden pressure jump corresponding to the low flow velocity. Accordingly, we described the stress wave propagation in total pressure and the flow velocity forms; we also visualized this process in the middle axial section (Figure 5). In this process, the dissipation effect was mainly due to the fluid viscosity and could not be neglected. We selected several points in the $z$ and $x$ directions on the bottom wall to figure out the proper areas for cell measuring (Figures 6–8). Only the $z$ shear stress changed slightly ( $6\%$ drop at every $15mm$ in this example), and the $x$ shear stress was so small that it could be neglected (in the interval of $0.25-0.75X$ ). These CFD results give us a high fault tolerance rate in repeating experiments.

Conclusions

In this paper, we focused on the application of dynamic loadings on single cells and revealed their mechanical response. Based on the Water-Hammer theories, we have developed a novel experimental method with a corresponding CFD model to help investigate cell mechanics under dynamic loadings. In this method, cells were normally cultured inside a microchannel. After impact, the stress wave generated applied a coupled compression-shear stress on an isolated living cell inside the microchannel. The results from an example model showed that the stress conditions could be easily controlled by controlling the velocity or length of the projectile. In addition, this method will allow researchers to adjust various design elements of their channels, such as the size, materials, etc., according to their lab’s environmental and actual needs, if the new design meets the relevant criteria presented in this study. This model offers repeatability, as a wide area is available for cell strain capturing, where cells suffer from nearly the same stress loadings.

Data availability statement

The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding authors.

Author contributions

YL, HY, NZ, and DX conceived and designed the research. DX, SW, and NZ derived the theoretical process and the CFD modeling process. DX, NZ, and PZ wrote the paper. All authors reviewed and edited the manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (12002285, 11722220, 61927810), Natural Science Basic Research Program of Shaanxi (2020JZ-11, 2020JQ-126), Fundamental Research Funds for the Central Universities (3102020smxy003), and the Program of Introducing Talents of Discipline to Universities (the 111 project) (No. BP0719007).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Ackroyd, J. (1967). On the laminar compressible boundary layer with stationary origin on a moving flat wall. Math. Proc. Camb. Phil. Soc. 63, 871–888. doi:10.1017/s0305004100041840

CrossRef Full Text | Google Scholar

Afzal, N. (1973). A higher order theory for compressible turbulent boundary layers at moderately large Reynolds number. J. Fluid Mech. 57, 1–25. doi:10.1017/s0022112073000996

CrossRef Full Text | Google Scholar

Alberts, B., Johnson, A., Lewis, J., Raff, M., Roberts, K., and Walters, P. (2002). Molecular biology of the cell. 4th edn. New York: Garland Science.

Google Scholar

Bao, G., and Suresh, S. (2003). Cell and molecular mechanics of biological materials. Nat. Mater 2, 715–725. doi:10.1038/nmat1001

PubMed Abstract | CrossRef Full Text | Google Scholar

Ben-Dor, G., Igra, O., and Elperin, T. (2000). Handbook of shock waves, three volume set. San Diego, CA: Elsevier.

Google Scholar

Bush, W. B., and Fendell, F. E. (1972). Asymptotic analysis of turbulent channel and boundary-layer flow. J. Fluid Mech. 56, 657–681. doi:10.1017/s0022112072002599

CrossRef Full Text | Google Scholar

Bush, W. B., and Fendell, F. E. (1973). Asymptotic analysis of turbulent channel flow for mean turbulent energy closures. Phys. Fluids 16, 1189–1197. doi:10.1063/1.1694497

CrossRef Full Text | Google Scholar

Davies, W., and Bernstein, L. (1969). Heat transfer and transition to turbulence in the shock-induced boundary layer on a semi-infinite flat plate. J. Fluid Mech. 36, 87–112. doi:10.1017/s0022112069001534

CrossRef Full Text | Google Scholar

Deshpande, V. S., Heaver, A., and Fleck, N. A. (2006). An underwater shock simulator. Proc. R. Soc. A 462, 1021–1041. doi:10.1098/rspa.2005.1604

CrossRef Full Text | Google Scholar

Desprat, N., Richert, A., Simeon, J., and Asnacios, A. (2005). Creep function of a single living cell. Biophysical J. 88, 2224–2233. doi:10.1529/biophysj.104.050278

PubMed Abstract | CrossRef Full Text | Google Scholar

Dudani, J. S., Gossett, D. R., Tse, T., and Di Carlo, D. (2013). Pinched-flow hydrodynamic stretching of single-cells. Lab. Chip 13, 3728–3734. doi:10.1039/c3lc50649e

PubMed Abstract | CrossRef Full Text | Google Scholar

Dymond, J., and Malhotra, R. (1988). The Tait equation: 100 years on. Int. J. Thermophys. 9, 941–951. doi:10.1007/bf01133262

CrossRef Full Text | Google Scholar

Evans, E., and Yeung, A. (1989). Apparent viscosity and cortical tension of blood granulocytes determined by micropipet aspiration. Biophysical J. 56, 151–160. doi:10.1016/s0006-3495(89)82660-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganadhiepan, G., Miramini, S., Patel, M., Mendis, P., and Zhang, L. (2019). Bone fracture healing under Ilizarov fixator: Influence of fixator configuration, fracture geometry, and loading. Int. J. Numer. Meth Biomed. Engng 35, (6), e3199. doi:10.1002/cnm.3199

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghidaoui, M. S., and Kolyshkin, A. A. (2001). Stability analysis of velocity profiles in water-hammer flows. J. Hydraul. Eng. 127, 499–512. doi:10.1061/(asce)0733-9429(2001)127:6(499)

CrossRef Full Text | Google Scholar

Guo, Z., Xu, K., and Wang, R. (2013). Discrete unified gas kinetic scheme for all Knudsen number flows: Low-speed isothermal case. Phys. Rev. E Stat. Nonlin Soft Matter Phys. 88, 033305. doi:10.1103/PhysRevE.88.033305

PubMed Abstract | CrossRef Full Text | Google Scholar

Hamed, D. (2020). Mechanical analysis of in vitro TBI models.

Google Scholar

Hao, Y., Cheng, S., Tanaka, Y., Hosokawa, Y., Yalikun, Y., and Li, M. (2020). Mechanical properties of single cells: Measurement methods and applications. Biotechnol. Adv. 45, 107648. doi:10.1016/j.biotechadv.2020.107648

PubMed Abstract | CrossRef Full Text | Google Scholar

Hoang, H., Galliero, G., Montel, F., and Bickert, J. (2015). Tait equation in the extended corresponding states framework: Application to liquids and liquid mixtures. Fluid Phase Equilibria 387, 5–11. doi:10.1016/j.fluid.2014.12.008

CrossRef Full Text | Google Scholar

Hutarew, G. (1973). Einführung in die Technische Hydraulik. Berlin, Heidelberg: Springer.

Google Scholar

Inaba, K., and Shepherd, J. E. (2010). Flexural waves in fluid-filled tubes subject to axial impact. J. Press. Vessel Technol. 132, 021302. doi:10.1115/1.4000510

CrossRef Full Text | Google Scholar

Joukowsky, N. (1900). Über den hydraulischen Stoss in Wasserleitungsröhren, Impériale des Sciences de St.Pétersbourg, Mémoires de lʾAcadémie.

Google Scholar

Kandlikar, S., Garimella, S., Li, D., Colin, S., and King, M. R. (2005). Heat transfer and fluid flow in minichannels and microchannels. New York, NY: Elsevier.

Google Scholar

Katti, D. R., Katti, K. S., Molla, S., and Kar, S. (2019). Biomechanics of cells as potential biomarkers for diseases: A new tool in mechanobiology. Encycl. Biomed. Eng. 13, 1–21. doi:10.1016/b978-0-12-801238-3.99938-0

CrossRef Full Text | Google Scholar

Korteweg, V. D. (1878). Ueber die Fortpflanzungsgeschwindigkeit des Schalles in elastischen Röhren. Ann. Phys. Chem. 241, 525–542. doi:10.1002/andp.18782411206

CrossRef Full Text | Google Scholar

Lennon, F. (1994). Shock wave propagation in water. Manchester, UK: Manchester Metropolitan University.

Google Scholar

Mandair, S. (2020). 1D and 3D Water-Hammer Models: The energetics of high friction pipe flow and hydropower load rejection. Canada: University of Toronto.

Google Scholar

Martins, N. M., Carriço, N. J., Ramos, D., and Covas, H. (2014). Velocity-distribution in pressurized pipe flow using CFD: Accuracy and mesh analysis. Comput. Fluids 105, 218–230. doi:10.1016/j.compfluid.2014.09.031

CrossRef Full Text | Google Scholar

Martins, N. M. C., Brunone, B., Meniconi, S., Ramos, H. M., and Covas, D. I. C. (2018). Efficient computational fluid dynamics model for transient laminar flow modeling: Pressure wave propagation and velocity profile changes. J. Fluids Eng. 140, 011102. doi:10.1115/1.4037504

CrossRef Full Text | Google Scholar

Martins, N. M., Soares, A. K., Ramos, H. M., and Covas, D. I. (2016). CFD modeling of transient flow in pressurized pipes. Comput. Fluids 126, 129–140. doi:10.1016/j.compfluid.2015.12.002

CrossRef Full Text | Google Scholar

Mellor, G. L. (1972). The large Reynolds number, asymptotic theory of turbulent boundary layers. Int. J. Eng. Sci. 10, 851–873. doi:10.1016/0020-7225(72)90055-9

CrossRef Full Text | Google Scholar

Melnik, R., and Grossman, B. (1974). Analysis of the interaction of a weak normal shock wave with a turbulent boundary layer. in 7th Fluid and PlasmaDynamics Conference, Palo Alto,CA,U.S.A., 598.doi:10.2514/6.1974-598

CrossRef Full Text | Google Scholar

Meyers, M. A. (1994). Dynamic behavior of materials. John Wiley & Sons.

Google Scholar

Mirels, H. (1955). Laminar boundary layer behind shock advancing into stationary fluid.

Google Scholar

Moeendarbary, E., and Harris, A. R. (2014). Cell mechanics: Principles, practices, and prospects. WIREs Mech. Dis. 6, 371–388. doi:10.1002/wsbm.1275

PubMed Abstract | CrossRef Full Text | Google Scholar

Nagayama, K., Mori, Y., Shimada, K., and Nakahara, M. (2002). Shock Hugoniot compression curve for water up to 1 GPa by using a compressed gas gun. J. Appl. Phys. 91, 476–482. doi:10.1063/1.1421630

CrossRef Full Text | Google Scholar

Nikpour, M., Nazemi, A., Dalir, A. H., Shoja, F., and Varjavand, P. (2014). Experimental and numerical simulation of water hammer. Arab. J. Sci. Eng. 39, 2669–2675. doi:10.1007/s13369-013-0942-1

CrossRef Full Text | Google Scholar

Patterson, L. H. (2020). The MicroHammer: Investigating cellular response to impact with a microfluidic mems device. Santa Barbara: University of California.

Google Scholar

Patterson, L. H., Walker, J. L., Rodriguez-Mesa, E., Shields, K., Foster, J. S., Valentine, M. T., et al. (2019). Investigating cellular response to impact with a microfluidic MEMS device. J. Microelectromechanical Syst. 29, 14–24.

Google Scholar

Pelham, R. J., and Wang, Y.-L. (1997). Cell locomotion and focal adhesions are regulated by substrate flexibility. Proc. Natl. Acad. Sci. U.S.A. 94, 13661–13665. doi:10.1073/pnas.94.25.13661

PubMed Abstract | CrossRef Full Text | Google Scholar

Pinnington, R. (1997). The axisymmetric wave transmission properties of pressurized flexible tubes. J. sound Vib. 204, 271–289. doi:10.1006/jsvi.1997.0938

CrossRef Full Text | Google Scholar

Rapp, B. E. (2016). Microfluidics: Modeling, mechanics and mathematics. Cambridge, US: William Andrew.

Google Scholar

Riasi, A., Nourbakhsh, A., and Raisee, M. (2009). Unsteady velocity profiles in laminar and turbulent water hammer flows. J. fluids Eng. 131, 121202. doi:10.1115/1.4000557

CrossRef Full Text | Google Scholar

Shepherd, J. E., and Inaba, K. (2009). Shock loading and failure of fluid-filled tubular structures. Dyn. Fail. Mater. Struct. 5, 153–190. doi:10.1007/978-1-4419-0446-1_6

CrossRef Full Text | Google Scholar

Smith, R. T. (1973). Weak shock wave propagation in liquid media.

Google Scholar

Thorley, A., and Guymer, C. (1976). Pressure surge propagation in thick-walled conduits of rectangular cross section.

Google Scholar

Trepat, X., Grabulosa, M., Puig, F., Maksym, G. N., Navajas, D., and Farré, R. (2004). Viscoelasticity of human alveolar epithelial cells subjected to stretch. Am. J. Physiology-Lung Cell. Mol. Physiology 287, L1025–L1034. doi:10.1152/ajplung.00077.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Urbanska, M., Muñoz, H. E., Shaw Bagnall, J. S., Otto, O., Manalis, S. R., Di Carlo, D., et al. (2020). A comparison of microfluidic methods for high-throughput cell deformability measurements. Nat. Methods 17, 587–593. doi:10.1038/s41592-020-0818-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Versteeg, H. K., and Malalasekera, W. (2007). An introduction to computational fluid dynamics: The finite method. London, UK: Pearson education.

Google Scholar

Walters, T. W., and Leishear, R. A. (2018). When the Joukowsky equation does not predict maximum water hammer pressures, in Pressure Vessels and Piping Conference, Prague, Czech Republic: American Society of Mechanical Engineers ASME. doi:10.1115/pvp2018-84050

CrossRef Full Text | Google Scholar

Yajnik, K. S. (1970). Asymptotic theory of turbulent shear flows. J. Fluid Mech. 42, 411–427. doi:10.1017/s0022112070001350

CrossRef Full Text | Google Scholar

Zhang, L., Gardiner, B. S., Smith, D. W., Pivonka, P., and Grodzinsky, A. (2008). A fully coupled poroelastic reactive-transport model of cartilage. Mol. Cell Biomech. 5 (2), 133–153. doi:10.3970/mcb.2008.005.133

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, L., Richardson, M., and Mendis, P. (2012). Role of chemical and mechanical stimuli in mediating bone fracture healing. Clin. Exp. Pharmacol. Physiol. 39 (8), 706–710. doi:10.1111/j.1440-1681.2011.05652.x

PubMed Abstract | CrossRef Full Text | Google Scholar