 Research
 Open Access
 Published:
A tracer liquid image velocimetry for multilayer radial flow in bioreactors
BioMedical Engineering OnLinevolume 14, Article number: 10 (2015)
Abstract
Background
This paper presents a Tracer Liquid Image Velocimetry (TLIV) for multilayer radial flow in bioreactors used for cells cultivation of tissue engineering. The goal of this approach is to use simple devices to get good measuring precision, specialized for the case in which the uniform level of fluid shear stress was required while fluid velocity varied smoothly.
Methods
Compared to the widely used Particles Image Velocimetry (PIV), this method adopted a bit of liquid as tracer, without the need of laser source. Subpixel positioning algorithm was used to overcome the adverse effects of the tracer liquid deformation. In addition, a neighborhood smoothing algorithm was used to restrict the measurement perturbation caused by diffusion. Experiments were carried out in a parallel plates flow chamber. And mathematical models of the flow chamber and Computational Fluid Dynamics (CFD) simulation were separately employed to validate the measurement precision of TLIV.
Results
The mean relative error between the simulated and measured data can be less than 2%, while in similar validations using PIV, the error was around 8.8%.
Conclusions
TLIV avoided the contradiction between the particles’ visibility and following performance with tested fluid, which is difficult to overcome in PIV. And TLIV is easier to popularize for its simple experimental condition and low cost.
Background
In cell cultivations for tissue engineering, the perfusion bioreactors with multilayer circular parallel plates scaffold were widely used. Ansede, Swift, Wang, Zhang, Lei et al. [15] all employed this type of perfusion bioreactors to implement their researches. However, relative few studies have addressed the uniformity of fluid shear stress among different layers. In fact, the magnitude of fluid shear stress significantly affects the proliferation and differentiation of cells [6,7]. And in many cell cultivations, the uniform level of fluid shear stress required was up to 0.001 Pa [811]. So it is necessary to determine the real magnitude of fluid shear stress in each layer by measuring the actual velocity distribution.
A widely used measurement method for fluid velocity is Particle Image Velocimetry (PIV) [12]. Jonathan, Bahar, Georgecet, Stephan, et al. [1316] used PIV to measure flows within bioreactors for cell culture purposes. PIV utilized advanced image acquisition, such as high speed cameras and laser light sources, and complicated computer processing technology to measure the fullfield flow. For it imaged in a crosssection illuminated with a sheet light source, it is unavoidable that some particles move into or shift out of the light sheet in measuring process. This change of particles resulted in the errors of algorithms based on correlation in PIV [17]. Moreover, the contradiction between the following performance and visibility of tracer particles is still difficult to overcome in PIV [18]. Reducing the size of particle improves the following performance with the fluid flow, but it weakens the visibility of particle [19]. The typical diameter of particle used in PIV is 3–300 μm. To balance this contradiction, Santiago; Chia et al. employed μPIV [20,21]. μPIV adopted 100–800 nm diameter microparticles and utilized microscope to take images. But such small particles significantly affect the accuracy of measurements [21]. And the small imaging region of microscope restricts the magnitude of velocity to only 5 mm/s currently [22]. Furthermore, the concentration of seeding particles maintained around 1% in PIV. It turned the original singlephase flow into solid–liquid twophase flow. In this situation, the consistency of the particles and flow field is of concern.
In the bioreactor with multilayer circular parallel plates scaffold, the flow field in every layer is radial flow. In cell cultivations, the goal of measurement is the distribution uniformity of fluid shear stress among different layers, but not the characterization of full flow field. This goal can be achieved by measuring the velocities of fluid among every layer at one vertical line. In a multilayer cylindrical radial flow field, when the total flow rate keeps constant, the position of the vertical line can be selected arbitrarily. In this research the focus is measurement precision, but not velocity gradient. Therefore, the measuring areas in which fluid velocity varied smoothly can be selected. For such measurement cases, a Tracer Liquid Image Velocimetry (hereinafter referred to as TLIV) was presented. The TLIV simply used general industrial digital cameras and natural light source, adopted liquid to replace particles as tracer. While for the tracer liquid is deformable, polydispersed and easy diffused, the accuracy of measurement is necessary to analyze. To eliminate the adverse effects of deformation, a subpixel positioning algorithm was utilized. And to overcome the measurement perturbations caused by the diffusion, a neighborhood smoothing algorithm was employed. For the lack of validated mathematical model of circular parallel plates bioreactors, the TLIV experiments were carried out in a rectangular parallel plates flow chamber. And to verify the precision, a mathematical model of the flow chamber and simulations using Computational Fluid Dynamics (CFD) [23] were performed, respectively.
Methods
Experimental setup
The experimental equipment for validating the TLIV was shown as Figure 1. Driven by a peristaltic pump, water flowed from a reservoir (unshown in Figure 1) into the tank in which a baffle was coupled to avoid affecting the flow field in parallel plates chamber. Then it flowed from inlet to outlet of the chamber driven by another peristaltic pump. Dye tracer was injected at the middle in the inlet through a syringe needle (unshown in Figure 1) with 0.5 mm aperture. The thickness and width of injected dye tracer should separately be kept within 1.5 mm and 5 mm to avoid adhering to the plates or lateral walls during it flowing through the chamber. A camera was placed right above or on the side of the parallel plates chamber to get images of the dye tracer flowing in the chamber.
Tracer liquid image velocimetry

i) Preprocessing
Before injecting tracer liquid, images of flow chamber filled with water were taken as the background images. Figure 2 A and B were the background images separately taken from top and side views. Let b(x, y) represent the intensity of background images, where (x, y) denote the coordinate of arbitrary pixel.
Then in all images, a rectangular region containing interesting objects was select as the processing windows, shown as the red frames in Figure 2.
For calibrating the scale in images, a graduate ruler was put into the middle layer in the flow chamber. Then the center lines of scales were extracted out. The number of pixels between every two center lines was the conversion coefficient k between pixel interval and actual size. Figure 2 CD illustrates the scale calibration.

ii) Subpixel positioning algorithm
Firstly, the acquired color images were transformed into gray images, notated as f(x, y), as shown in Figure 2 F.
To enhance the interesting objects, a subtraction operation was carried out by using the form
where b(x, y) represents the intensity of background images.
Then a median filter [24] was used to reduce noise. Afterward, to extract the objects, an intensity thresholding algorithm was implemented by the expression
where T is an intensity threshold and f ' (x, y) is the output images of the median filter.
To determine the boundary, B, of interesting object, an edge detection algorithm with Canny edge detector [24] was carried out. Afterward, the coordinate of tracer mass center can be calculated by using the following expression
where (x _{ i }, y _{ i }) denote the coordinate of arbitrary pixel i, for all i ∈ B, and f _{ i }(x _{ i }, y _{ i }) denote the intensity of pixel i. The Eq. (3) is of subpixel algorithm that the coordinate obtained can be less than 1 pixel.
Figure 2 G15 illustrates the procedures of the subpixel positioning algorithm.

iii) Postprocessing
Let set A represent the extracted coordinates of object positions from a set of continual images. Then the instantaneous displacement speed of objects between any two continual image i and image i1 can be calculated by
where a _{ i }, a _{ i − 1} ∈ A and t _{ i }, t _{ i − 1} separately denote the photographing time (in unit of millisecond) of images i and i1.
For reducing the perturbations of the displacement, a neighborhood smoothing [24] algorithm was implemented by using following expression
where r denote the neighborhood range of object position i.
Then the moving velocity of object at a _{ i } can be obtained by
where 1/k is the physical distance ( in unit of millimeter) between two neighboring pixels, obtained in preprocessing.
Figure 2 H12 illustrates the procedures of postprocessing.
According to the kinematic law of mass center in theoretical mechanics [25], in mass system, the movement law of mass center is completely equal with the movement of one mass point. The acceleration of mass center depends on the magnitude and direction of principal vector of the external force system, while is independent of internal force of mass system and the acting position of external force. Though the tracer liquid cannot be rigid, and it will be deformed and diffused in the flow field. If the diffusion exhibits symmetry statistically, the effect of diffusion will be counterbalanced. The time between two continuous images of tracer was less than 60 ms (for actual frame rate ≥ 17 fps). Consequently, for incompressible fluid, if the deformation of tracer liquid profile can be ignored during such short time, the tracer can be considered as rigid body. Then the theory mentioned in reference [25] can be used for tracer liquid. It means that the moving velocity of mass center can represent the tracer velocity. For verifying these assumptions, the validating experiments were carried out in this paper.
Experimental validations and operating conditions
The parallel plates flow chamber shown as Figure 1 A was 3 mm high, 38 mm wide and 55 mm long. Tested fluid was liquidwater. The rotational speed of the peristaltic pump was 200 rpm and the flow out rate was 57 ml/min. The image processing windows A and B shown in Figure 2 AB separately were 370 pixels and 50 pixels high. Two windows both were 490 pixels wide, and the distances between each window and the outlet boundary of flow chamber both were 2 mm.
Tracer liquid was waterdiluted black ink (Hero $$ R $$ ) with a density of 0.999 g/cm^{3}. The injecting position of tracer located at the horizontal middle layer of inlet. The injecting direction was normal to flowing direction of tested fluid.
Imaging device was industrial digital camera (CatchBEST $$ R $$ MU3C500C) with USB 3.0 interface and frame cache. Imaging background was white baffle illuminated by natural light source. The frame rate was 15 – 20 fps @ 720 × 518 pixels. Images were processed by software MATLAB $$ R $$ .
Photographing positions were separately located at top and side of parallel plates flow chamber.
Comparing validations
To verify the TLIV experiments, mathematical models of the parallel plates flow chamber were utilized. The velocity model is [26]
where W, H, Q _{ V }, y separately represent the width, height, flow rate of flow chamber, and the vertical coordinate. The original point located at the center of inlet. When y = 0 mm, the velocity gets its maximum value
where ū is mean velocity
The fluid shear stress on parallel plates is
where μ is viscous coefficient (this paper μ = 0.001003 kg/ms).
Moreover, Computational Fluid Dynamics was used to simulate the flow field. The size of simulated geometry was the same as the parallel plates flow chamber. The geometry was meshed in software ICEM (ANSYS $$ R $$ ) with nonstructure grid. The flow field was calculated in software Fluent (ANSYS $$ R $$ ) with inlet velocity of 0.0083 m/s, steady laminar model ignoring energy transform.
Results
Separately according to Eq. (9) and Eq. (8), the mean and maximum velocities separately were 8.3 mm/s and 12.5 mm/s at the flow rate of 9.5E7 m^{3}/s. On the base of Eq. (7), for y = −0.5 mm and −0.6 mm, the velocities were 11.11 mm/s and 10.5 mm/s, respectively.
By CFD, The flow rate was 9.4996E7 m^{3}/s. And for y = −0.5 mm and −0.6 mm, the distributions of velocity vector were shown as Figure 3.
By the preprocessing shown as Figure 2 CD, the value of coefficient k in Eq. (6) was 15.5 pixel/mm.
The illustrations of TLIV on top view were shown in Figure 4. Its velocity distribution at 5 observing points were shown in Table 1. In Figure 4, M15, U15 and D15 were corresponding to Cases 2, 1, and 3 in Table 1, respectively.
The illustrations of TLIV at side view were shown as Figure 5. Its velocity distributions were shown in Figure 6. For comparing, the simulated data was also plot in Figure 6.
For 3 trial cases of TLIV at side view, the velocities at 5 observing points were shown in Table 2. Figure 5 was corresponding to the Case 1 in Table 2.
The variable C _{y} shown in Table 2 can be transformed to y by following expression
where 3 represented the height (in unit of mm) of the chamber, multiplying by 15.5 to transform its unit of mm into pixel. Then Eq. (11) indicated the calculation of the vertical distance (in unit of mm) between the position descripted by C _{y} and the horizontal middle layer in the chamber.
Substituting the values of y into Eq. (7), the modelbased velocity, V _{C}, can be calculated. The errors between V shown in Table 2 and V _{C} can be obtained and were shown as Table 3.
The variable C _{x} can be transformed to L shown in Figure 3 by expression
where 55 was the length (in unit of mm) of the chamber. 490 was the length (in unit of pixel) of processing window. And 2 was the distance (in unit of mm) between the left side of processing window and the left side of the chamber. Then Eq. (12) indicated the calculation of the horizontal distance (in unit of mm) between the position descripted by C _{x} and the inlet of the chamber.
Based on Eqs. (11) and (12), the coordinate (C _{x}, C _{y}) in image processing window can be transformed to the coordinate of (L, y) shown in Figure 3 where y indicated the vertical height in flow chamber. Then the simulated velocity of flow field corresponding to (C _{x}, C _{y}) can be collected from simulated results through indexing (L, y).
As a result, combining the values of y, L, and W = 0, the simulationbased velocity, V _{CFD}, can be collected. Then the errors between V and V _{CFD} can be obtained and were shown as Table 4.
Discussion
According to Eq. (7), the velocity distribution is independent of fluid viscosity. So in this paper the water was used for validation.
In Figure 3, the fluid near outlet flowed from two sides toward middle. This phenomenon was caused by the stagnant flow in two corners. The pathlines of mass center shown in Figure 4 U5 and D5 exhibited the same result. Thus the top view of TLIV reflected the stagnant flow existing in the flow chamber.
For the stagnant, the actual width of flow field narrowed as fluid flowed to outlet. When total flow rate keeps constant, velocity increases as the width of flow field decreases. Figure 3 and Tables 1 and 2 indicated that the simulated and measured results both exhibited this trend. But Eq. (7) is based on the assumption that the length of flow chamber is infinite. It can’t reflect horizontal velocity variation.
In Table 1, some values deviated the mean velocity more than 1 mm/s. These conspicuous errors were resulted by the sharp change of background illuminant shown as Figure 2 A. While in Figure 2 B the illumination of background image was smoother. So in Table 2, at each observing point, the errors between every velocity and the mean velocity were less than 0.4 mm/s. At side views of the flow chamber, for the light diffused through a longer path, sharp illuminant changes were avoided. Thus Good performance of TLIV can be achieved under natural light source. The following discussions were mainly based on side views.
In all experiments, tracer was injected at horizontal middle layer of inlet. For the density of tracer is close to which of water, the vertical height of tracer can keep steady in flow chamber. The values of C _{y} in Table 2 meet this characteristic. Its variation was less than 3 pixels.
In Figures 4 and 5, the deformations of tracer liquid were obvious. Such as from M1 to M5 in Figure 4, the profile of tracer varied toward lengthening. From Eq. (7) it can be seen that the velocities of flow in different position within the chamber exist discrepancy, which resulted in the obvious deformations of tracer. But from their mass centers extracted in Tables 1 and 2, it can be seen that the subpixel positioning algorithm can eliminate the adverse effects of deformations.
A noticeable effect on object positioning is liquid diffusion, shown as the black curve in Figure 6. The red curve in Figure 6 demonstrated that the neighborhood smoothing algorithm can effectively restrict such adverse effects. However, the smoothing process simultaneously filtered the sharp changing components of velocity. As shown in Figure 6, at the positions of column = 170, 370, and 400 pixels, though the simulated values changed obviously by the variation of C _{y}, the measured values didn’t follow these abrupt changes. This discrepancy was resulted by the influence of smoothing algorithm. For such sharply varied velocity, the maximum relative error between measured and simulated values was up to 16.4%, shown as positions 4–5 in Table 4. But as mentioned before, the suitable positions for measurement using TLIV can be selected in the region where the velocity keeps smooth, such as positions 1–3 shown in Table 2. In these areas, the mean relative error range between measured and simulated data was [1.5%,1.8%]. While in similar validations by PIV, the error was around 8.8% [27,28].
According to Eq. (7) and Eq. (10), the fluid shear stress is proportion to velocity. Therefore, the distribution of fluid shear stress can be determined by the measurement of velocity.
For the dilution effects of water on the tracer liquid, the TLIV is not applicable to the situations that need to trace over a long period of time. All tests in this paper were no more than 5 seconds. In next research, we will find a tracer liquid that can’t be dissolved by water.
Though there is discrepancy between the flow field in rectangular flow chamber and that in circular parallel plates bioreactors, the main purpose of this paper is to validate the feasibility and measuring precision of TLIV. In future study, we will establish mathematical model for circular parallel plates bioreactors, and experimentally verify it with the TLIV.
Conclusions
In cell cultivations, it is necessary to measure the velocity distribution of fluid in bioreactors. For a cylindrical multilayer radial flow field, we can select the measuring regions in which velocities varied smoothly to determine the distribution uniformity. To get good measuring precision by simple devices, this paper provides a velocimetry. It used tracer liquid to replace particles, and employed digital image processing technology to solve the problems caused by the new tracer. Its measuring precision sufficiently met the requirement of application.
Abbreviations
 CFD:

Computational Fluid Dynamics
 PIV:

Particles Image Velocimetry
 TLIV:

Tracer Liquid Image Velocimetry
References
 1.
Ansede JH, Smith WR, Perry CH, Brouwer KR. An in vitro assay to assess transporterbased cholestatic hepatotoxicity using sandwichcultured rat hepatocytes. Drug Metab Dispos. 2010;38(2):276–80.
 2.
Lei X, Talha A, Zhang SF, Tuo XY. Hepatocyte function within a stacked double sandwich culture plate cylindrical bioreactor for bioartificial liver system. Biomater. 2012;33:7925–32.
 3.
Swift B, Pfeifer ND, Brouwer KL. Sandwichcultured hepatocytes: an in vitro model to evaluate hepatobiliary transporterbased drug interactions and hepatotoxicity. Drug Metab Rev. 2010;42(3):446–71.
 4.
Wang DL. Stem cells tissue engineering technology: basis theory and clinical application. China: Scientific Publishers; 2011. p. 669.
 5.
Zhang Y, Shi XL, Han B, Gu JY, Chu XH, Xiao JQ, et al. Immuno safety assessment of a novel multilayer flatplate bioartificial liver based on nanofiber scaffolds. Chin J Gen Surg. 2012;21:41–7.
 6.
Saini S, Wick TM. Concentric cylinder bioreactor for production of tissue engineered cartilgae: effect of seeding density and hydrodynamic loading on construct development. Biotechnol Prog. 2003;19(2):510–21.
 7.
Sebastian DB, Silvia T, Sezin EO, Toon L, Hans VO, Daniel B, et al. Bimodular flow characterization in tissue engineering scaffolds using computational fluid dynamics and particle imaging velocimetry. Tissue Eng: Part C. 2010;16(6):1553–64.
 8.
Bancroft GN, Sikavitsas VI, Dolder J. Fluid flow increases mineralized matrix deposition in 3D perfusion culture of marrow stromal osteoblasts in a dosedependent manner. Proc Natl Acad Sci U S A. 2002;99(20):12600–5.
 9.
Cherry RS, Papoutsakis ET. Physical mechanisms of cell damage in microcarrier cell culture bioreactors. Biotech Bioeng. 1988;32(8):1001–14.
 10.
Godwin TJ, Prewett TL, Wolf DA. Reduced shear stress: a major component in the ability of mammalian tisues to form threedimensional assemblies in stimulated microgravity. J Cel Biochem. 1993;51:301–11.
 11.
Porter B, Zauel R, Stockman H, Guldberg R, Fyhrie D. 3D computational modeling of media flow through scaffolds in a perfusion bioreactor. J Biomech. 2005;38(3):543–9.
 12.
Willert CE, Gharib M. Digital particle image velocimetry. Exp Fluids. 1991;10:181–93.
 13.
Bahar B, Gilda A. Location of scaffolds in bioreactors modulates the hydrodynamic environment experienced by engineered tissues. Biotechnol Bioeng. 2007;98(1):282–94.
 14.
Georgec E, Sarahc V, Stephanusg B, William JF, Krishnanb C, Davida V, et al. A novel flexstretchflow bioreactor for the study of engineered heart valve tissue mechanobiology. Ann Biomed Eng. 2008;36(5):700–12.
 15.
Jonathan D, John S, Kerry H. A fluid dynamics approach to bioreactor design for cell and tissue culture. Biotechnol Bioeng. 2006;94(6):1196–208.
 16.
Stephan CK, Valentin J, Carmen S, Dieter E, Silke B, Christian VB, et al. Fluid flow and cell proliferation of mesenchymal adiposederived stem cells in smallscale, stirred, singleuse bioreactors. Chem Ing Tech. 2013;85(1–2):95–102.
 17.
Hall JF, Barigou M, Simmons MJH, Stitt EH. A PIV study of hydrodynamics in gas–liquid high throughput experimentation (HTE) reactors with eccentric impeller configurations. Chem Eng Sci. 2005;60:6403–13.
 18.
Gao DR, Wang YQ, Shen GX. DPIV Technique and its application in flow field measurement. Hyd Pneum Seals. 2001;89:30–3.
 19.
Zhao Y. The study of tracer’s following behaviors in particles image velocimetry measurements. China: Dalian University of Technology; 2004. p. 51.
 20.
Chia ML, Abram V, Gary B, Timothy W. Flow bioreactor design for quantitative measurements over endothelial cells using microparticle image velocimetry. Rev Sci Instrum. 2013;84:045109.
 21.
Santiago JG, Wereley ST, Meinhart CD, Beebe DJ, Adrian RJ. A particle image velocimetry system for microfluidics. Exp Fluids. 1998;25:316–9.
 22.
Lin TJ, Chen PC. Studies on hydrodynamics of an internalloop airlift reactor in gas entrainment regime by particle image analyzer. Chem Eng J. 2005;108:69–79.
 23.
Shouliang Q, Zhenghua L, Yong Y, Han JW, Yan K. Computational fluid dynamics simulation of airflow in the trachea and main bronchi for the subjects with left pulmonary artery sling. BioMed Eng OnLine. 2014;13:85.
 24.
Gonzalez RC, Woods RE. Digital image processing. China: Publishing House of Electronics Industry; 2010. p. 241–53.
 25.
Liu YZ, Yang HX, Zhu BH. Theoretical mechanics. China: Higher education press; 2001. p. 206–7.
 26.
Liu ZR, Jiang WY, Huo Y. The analysis of the steady flow in the flow chamber. J Hydrod. 1997;2:37–45.
 27.
Wang X, Ding J, Guo WQ, Ren NQ. A hydrodynamics–reaction kinetics coupled model for evaluating bioreactors derived from CFD simulation. Biores Tech. 2010;101:9749–57.
 28.
Niu LL, Qian M, Wan K, Yu WT, Jin QF, Ling T, et al. Ultrasonic particle image velocimetry for improved flow gradient imaging: algorithms, methodology and validation. Phys Med Biol. 2010;55:2103–20.
Acknowledgements
The authors appreciate the support from the National Natural Science Foundation of China (No. 51205421), the Postdoctoral Science Foundation of China ( No. 2012M521647), the Science and Technology Program of Guangzhou (China [2013]164), and the Key Laboratory of Sensing Technology and Biomedical Instruments of Guangdong province (China 2011A060901013).
Author information
Additional information
Competing interests
The authors declare that they have no competing interests.
Authors’ contributions
YG conceived, validated the Tracer Liquid Image Velocimetry and drafted the manuscript. JL simulated the flow fields. YL supervised the project and research group. JY drew the figures and tables. All authors read and approved the final manuscript.
Rights and permissions
About this article
Received
Accepted
Published
DOI
Keywords
 Bioreactor
 Fluid shear stress
 Image velocimetry
 Computational Fluid Dynamics
 Parallel plates flow chamber