Team:SUSTech Shenzhen/Model

Model

Overview

Modeling is usually used to make sense of the experimental discovery in traditional biological studies. In this synthetic biology project, we believe that carefully carried out modeling will be critical for the experiment design and data analysis at different stages of the project. We hope to demonstrate that the modeling is especially helpful to finalize the microfluidics chip design, and determine the absolute shear force in the experiments, as well as characterization and optimization of acoustic stimulation experiments.

Simplified fluidic mechanic theory is sufficient to guide the geometry design of shear force screening microfluidics chip. Together with direct experimental measurement of flow speed, numerical simulation of the flow in the microfluidics chip using Finite Element Analysis (FEA) enable us to derive actual shear force in the chip.

On the other hand, the acoustic stimulation devices are very different from each other. Our empirical is not sufficient to help select and data analysis. We use Finite Element Analysis (FEA) to model the force field generated by the acoustic stimulators, with the help of acoustometer measurements.

Force Field

To quantitatively measure the response of MS channel on the cell membrane, we designed our experiment manipulated on microfluidic chip. CHO-K1 cells were cultured adherent on the bottom of PDMS channels. Stable and flexible force field could be formed surrounding cell, which directly controlled by the pumped-inflow rate. Each time when we applied a constant pumped inflow, cells in 3 different observation channels could receive corresponding small, middle, large, 3 level of force magnitude(1: 9: 81). By changing the pumped flow rate, we could measure MS channel response under a series of forces with different magnitude order.

Fig. 1 Microfluidic Channels

Mathematical Analysis

During experiment, when fluid(culture media x flowed through the channels, the shear force was applied on a thin layer close to wall. The culture media can be treated as a Newtonian flow with a constant viscosity μ,0.012dyn·s/cm². The shear stress that medium generated on the wall is proportional to strain rate. As a result, we could relate the magnitude of stress to the velocity gradient along the transversal surface of each channel.

Due to the extremely small scale of microfluidic channel (only 0.285mm of width and 0.090mm of height), fluid in it exists as the Laminar flow. Turbulent flow is not possible with the physiology range of flow speed, even at the corners. The head loss of the pumped inflow is also negligible in our experiment. The modeling of flow between 2 parallel plates applies to such stable state, in which the magnitude of shear stress on the bottom is proportional to the maximum velocity ( shown as the longitudinal central line in Fig. 2).

As the shear stress is proved proportional to the flow rate in the channels, we designed the chip to achieve that each of 3 channels in observation area has the same transversal surface but the gradient (81: 9: 1 from up to bottom) flow rate, also the average speed, as shown in Fig. 3.

Fig. 2 Speed distribution inside channels

Newtonian fluid flows between two wide, parallel plates Flow driven by pressure difference Parabolic velocity profile given by u = V[1-(y/h)^2] where V is the maximum velocity (along channel centerline y=0) Along bottom wall, y=-h, shear stress \tau=\mu \frac{du}{dy}=\mu(-2V\frac{y}{h^2})=\frac{2\mu V}{h} For a constant flow rate within the channels, the average speed remains same \bar V=\frac{Q}{A}=\frac{\int_A udA}{A}=\frac{1}{2}V

Fig. 3 The geometry of the microfluidic channels. The lengths of straight channels between AB and AC (L_AB, L_AC) are 14.3mm and 12.0mm, and the curved channels between AB and BC(R_AB, R_BC) are 108.0 mm and 105.0mm.

Based on the integral form of the continuity equation for an incompressible fluid, at any given time, the sum of flow volume in each branch is equal to the total volume pumped in the whole chamber. According to the Kirchhoff's law, the accumulation of flow speed at the joint point must be zero. It is easy to get the equation that\begin{matrix} L_{AC}V_1=L_{AB}V_2+R_{AB} V_3 \\ V_2:V_3=R_{AB}:L_{AB} \end{matrix}

Consequently, we could reach the conclusion that the maximum velocity, average flow rate and the shear force on the bottom of each channel is at a ratio of 81 : 9 : 1. To verify the modeling result, we made calibration in the real microfluidics chip, and the testing result is very close to designing value (about 75 : 8 : 1 of the maximum velocity ratio in each channel) and the simulation result (80 : 8 : 1).

Experimental Result Analysis

In the whole experiment, we tested one mechanosensitive channel, Piezo1, with different magnitude orders of shear stress.

Comparing the fluorescence intensity amplitudes corresponded to the shear stress ranged from 0 Pa to 5 Pa, we found a best response at around 0.01-0.1 Pa. No response was observed when the shear stress is larger than 2 Pa.

Fig. Comparison between the R-GECO cell response(a) and Piezo1+R-GECO cell response(b) to the 15kHz Buzzer’ s simulation.

Fig. 5 Fluorescence change of R-GECO+Piezo1 cell when stimulate by 20kHz sound emitted by a buzzer. Sound was started from 5s to 15s. 5 cells were selected and their fluorescence intensity changed as shown in the figure.

Discussion Details

Calibration

Rainbow beads (SPHEROTM Rainbow Calibration Particles, diameter: 6um) were dissolved in the cell culture medium. A series of pumped flow rate (5ul/min, 15ul/min, and 45ul/min) were applied to generate a steady fluid flow. The exposure time was set to 100ms (one representative image is shown in Fig. 6).

Fig. 6 Bead trace in the fastest microfluidics channel with The Flow Rate of 5ul/min

By using our customerized MATLAB program, the ratio of maximum flow velocities were calculated as 13 : 100 : 978. The absolute speed with pump speed of 45 um/min is shown in Fig. 7.

Fig. 7 Maximum Velocities in Each Channel with The Flow Rate of 45ul/min

Code & Details

Simulation

To further validate our design, we also performed the simulation using the Finite Element Method (FEM) to confirm that the field in microfluidic channel satisfies our expectation. FEM is a method for modeling and simulating physics-based problems. It divides a large problem into smaller elements, these smaller elements will be solved separately and then be assembled backed to the original model. A collection of smaller elements is called a “mesh”, the “mesh” is constructed by numerous elements, the smaller the size of the element is, the finer the “mesh” is, and the more accurate will the result be. The final geometry set up can be seen from below (Fig. 8).

Fig. 8 Final geometry setup

Under different pumped inflow rate(5μL/min, 15μL/min and 45μL/min), the velocity magnitude distribution was calculated. The result of maximum velocity and shear force applied on the bottom wall matched with the mathematics perfectly(around 80: 10: 1 in all the conditions, shown in the table below).

Channel Flow Rate	fast	middle	slow
5μL/min	5414μm/s	549μm/s	68μm/s
15μL/min	16242μm/s	1646μm/s	204μm/s
45μL/min	48724μm/s	4938μm/s	612μm/s

Table 1 Simulation results of central line speed of each channel

Results of mathematical demonstration (1: 9: 81), calibration (13: 100: 978) and simulation (1: 10: 80) are consistent with each other. This shows that the model is feasible for estimating the shear force that generated by microfluidics on cell membrane. In order to obtain the exact shear force, we inputted the experimentally used flow rate as a boundary condition (6.25μL/min,50μL/min) in our model. The shear stress was also calculated.

Channel Flow Rate	fast	middle	slow
6.25μL/min	0.269Pa	0.027Pa	0.003Pa
50μL/min	2.15Pa	0.22Pa	0.03Pa

Table 2 Simulation results of central line shear stress

Fig.9 Shear stress distribution under input flow rate 50μL/min (log scale)

Detailed Protocol

Acoustic Field

Sonic and ultrasonic waves could generate a pressure field to the cells adhering to the bottom of the cell dish. In our study of the sensitivity of mechanosensitive (MS) channel, we really found a most efficient sound frequency using buzzer and 108kHz ultrasonic transducer.

Fig. 10 B

Fig. 10 A

Fig. 10 (A) R-GECO+Piezo1 cell fluorescence intensity increased greatly when stimulated by 15K and 20KHz sound. While no obvious fluorescence intensity changed when GECO cells (without MS channel) were stimulated (B).

However, to quantitively define the relationship of channel opening ( considering the influence of energy or force difference), it was important to know the pressure field distribution on the cell layer under certain conditions. Mathematical analysis and finite element analysis (FEA) were used to define the modeling.

Audible Frequency Result Analysis

To discover the response of MS channels to audible vibration, buzzers, balanced amatures and speakers were tried as sonic stimulators applied to the cells with R-GECO and Piezo1 channel, but only buzzers induced obvious cell response significantly in high frequency conditions (15kHz and 20kHz) when they are attached to the bottom of a cell dish each time.

To explain that, a modal analysis method could be used when a buzzer was attached to the bottom of a cell dish, the shear stress induced by the asynchronism between cell culture medium and cell dish buzzer-forced oscillation is the driving factor of the channel response. It is interesting that the shear force magnitude in this condition is similar to the force magnitude in microfluidics channels.

Fig. 12 Buzzer attached on the bottom of cell dish

Fig. 11 Design pattern of the buzzer


Considering that the cell culture medium is a viscous and incompressible fluid with μ= 0.012dyn·s/cm2(0.0012Pa·s), it can not vibrate uniformly with the buzzer. The shear stress due to the differential distribution of fluidic velocity surrounding cell can apply a strain \gamma on the cell membrane.

\dot \gamma=\lim_{\delta t \rightarrow 0}\frac{\delta \theta}{\delta t}=\frac{\partial \theta}{\partial t}=\frac{\partial u}{\partial v}

For sound waves, I \propto p_0^2, p₀ is the acoustic pressure and I \propto s_0^2, s₀ is the displacement amplitude.

Thus we have I = \frac{p_0^2}{2\rho v},\Delta p_m=(v\rho \omega)s_m, which indicates that the intensity of sound waves follow an inverse square law.

Shear stress in the fluid is given by \tau = \mu \frac{du}{dy},u=s_m\cdot f

To define the amplitude of dish vibration, we measured the sound intensity by using Brüel & Kjær 1/8-Inch Pressure Field Microphone with Type 2669 Preamplifier (4138-C-006), which is calibrated from 6 Hz to 140 kHz.
 Sound speed in the air v=340 m/s, air density is 1.29 kg/m³, angular velocity is 2π*15000Hz forced by the buzzer.

Fig. 13 Measure the sound intensity by acoustometer

Acoustic pressure \Delta p_m = 4\text{Pa} is applied by the buzzer attached on the bottom of dish. 
From the formula above, the amplitude of dish vibration is around 9.68*10-7m. The shear stress applied on the dish bottom could also be calculated and the result is 0.01 Pa, which matches accurately to the shear stress generated in the microfluidic channel when the flow rate is 50μL/min in the slow channel.

Ultrasound Experiment Analysis

First of all, a function was studied to give pressure on the surface of the ultrasonic transducer-water interface knowing the input voltage. Based on the interface pressure and physical properties of device, we constructed FEA models to simulate the ultrasonic radiant pressure field of the cell layer. Pressure field was studied in conditions of diverse power and distance between transducer interface and cell layer.

Fig. 14 R-GECO+Piezo1 cells with DOX added in culture medium were stimulated by 108kHz ultrasound from 10s to 30s. Four regions were selected. Two R-GECO+Piezo1 cell region were chosen one from the right side of the video and another from left side. Two backgrounds were chosen the same way.

Fig. 15 Equivalent Schematics of the Transducer

The equivalent circuit of ultrasonic transducer was used to calculate the ultrasonic radiant power of an ultrasonic transducer. We assume the ultrasonic energy was evenly distributed on the transducer interface, that is:

I\left( x,y,z = 0 \right) = \left\{ \begin{matrix} I_{0},x^{2} + y^{2} \leq R^{2} \\ 0,x^{2} + y^{2} \lt R^{2} \end{matrix} \right.\

Where I₀ is the absolute ultrasound intensity, R is the radio of transducer interface

According to power intensity and pressure relation:

I_{0} = \frac{{P_{0}}^{2}}{2c\rho}

In which P_0 , c, \rho are respectively absolute pressure, acoustic velocity and water density. We managed to calculate the pressure on the transducer-culture medium interface.

Then, we implemented FEM to do a series of simulations to find the absolute pressure filed distribution in different conditions. 108kHz and 1.1MHz transducer device were modeled as is shown in Fig. 16

Fig. 16 B

Fig. 16 A

Fig. 16 D

Fig. 16 C

Fig. 16 (A) Pressure distribution of 108KHz ultrasonic transducer, the transducer was 1.0mm far above the cell layer. (B) Pressure distribution of 1.1MHz ultrasonic transducer. 3D pressure distribution for 108KHz condition (C) and 1.1MHz condition (D).

The influence of distance error on pressure distribution pattern was also studied, which indicated that the 108kHz transducer was very robust to distance error, while the 1.1MHz acted as a chaos when little distance error occurs.

In all, the 108kHz ultrasonic transducer was recommended as ultrasonic stimulator for quantitative measure cell response based on modeling and experimental results.

Details

References

Made by from the iGEM team SUSTech_Shenzhen.

Licensed under CC BY 4.0.

• ©2016 SUSTech_Shenzhen