EVALUATION OF BLOOD PLASMA FLOW AROUND 2D ERYHTROCYTES WITH THE USE OF COMPUTATIONAL FLUID DYNAMICS

INTRODUCTION Hemodynamic forces, such as Wall Shear Stress, are known to be one of the factors behind atherosclerotic plaque formation in blood vessels. This plaque formation can lead to vascular pathologies, such as aneurysms and stenosis (BALOSSINO et al., 2009; CEBRAL et al., 2005; DONG; INTHAVONG; TU, 2013). Given the importance of understanding the hemodynamics inside blood vessels, numerical simulations such as those using the CFD (Computational Fluid Dynamics) technique have been applied, in order to aid physicians with further analysis of the blood flow. These simulations can be used alongside with common medical imaging technologies, such as MRI (Magnetic resonance Imaging) and CT-scans (Computed Tomography) (LI et al., 2018). Although numerical studies try to simulate cases as close as possible to the real physical and biological behavior, many studies tend to adopt simplifications. Some usual considerations are the adoption of a Newtonian model for the single-phase blood rheology or the use of generic geometries representing the vessels. Nevertheless, the use of simplifications can lead to reasonable results (KARIMI et al., 2014; MORALES et al., 2013), which can be used as basis for further studies. Despite achieving acceptable results, blood is a suspension-type fluid, composed of blood cells, erythrocytes, leucocytes and thrombocytes, suspended in plasma (JUNG; HASSANEIN, 2008; RAZAVI; SHIRANI; SADEGHI, 2011). Whilst blood cell presence in blood flow is a known information, simulations tend to consider blood as a single-phase flow. Therefore, computational hemodynamics tend to focus on the whole blood as well as pathologies, rather than the physics behind the cells. However, blood cells can account to more than half of the blood volume, depending on the patient (JUNG; HASSANEIN, 2008; OU et al., 2016). Hence, the present work aimed to study the behavior of plasma, a Newtonian fluid, flowing around a single erythrocyte, as well as a group of erythrocytes contained in a 2D domain. From the simulations, parameters such as velocity magnitude and Reynolds number were analyzed, in order to provide more information that will be used in future studies by the research group.


INTRODUCTION
Hemodynamic forces, such as Wall Shear Stress, are known to be one of the factors behind atherosclerotic plaque formation in blood vessels. This plaque formation can lead to vascular pathologies, such as aneurysms and stenosis (BALOSSINO et al., 2009;CEBRAL et al., 2005;DONG;INTHAVONG;TU, 2013). Given the importance of understanding the hemodynamics inside blood vessels, numerical simulations such as those using the CFD (Computational Fluid Dynamics) technique have been applied, in order to aid physicians with further analysis of the blood flow. These simulations can be used alongside with common medical imaging technologies, such as MRI (Magnetic resonance Imaging) and CT-scans (Computed Tomography) (LI et al., 2018). Although numerical studies try to simulate cases as close as possible to the real physical and biological behavior, many studies tend to adopt simplifications. Some usual considerations are the adoption of a Newtonian model for the single-phase blood rheology or the use of generic geometries representing the vessels. Nevertheless, the use of simplifications can lead to reasonable results (KARIMI et al., 2014;MORALES et al., 2013), which can be used as basis for further studies. Despite achieving acceptable results, blood is a suspension-type fluid, composed of blood cells, erythrocytes, leucocytes and thrombocytes, suspended in plasma (JUNG; HASSANEIN, 2008;RAZAVI;SHIRANI;SADEGHI, 2011). Whilst blood cell presence in blood flow is a known information, simulations tend to consider blood as a single-phase flow. Therefore, computational hemodynamics tend to focus on the whole blood as well as pathologies, rather than the physics behind the cells. However, blood cells can account to more than half of the blood volume, depending on the patient (JUNG; HASSANEIN, 2008;OU et al., 2016). Hence, the present work aimed to study the behavior of plasma, a Newtonian fluid, flowing around a single erythrocyte, as well as a group of erythrocytes contained in a 2D domain. From the simulations, parameters such as velocity magnitude and Reynolds number were analyzed, in order to provide more information that will be used in future studies by the research group.

MATERIALS AND METHODS
In order to carry out simulations, three generic geometries were developed in the pre-processing step. The first two geometries had a single erythrocyte, where the differences between both are the orientation and shape of the cell. In the first one, the erythrocyte is oriented parallel to the flow, whereas the second one is oriented perpendicular to the flow. The third geometry is composed of twenty-one erythrocytes oriented parallel to the flow. Data for superficial area, as well as erythrocyte volume was taken from the literature (JAY, 1975). With the geometrical data, 3D erythrocytes were made, in order to assure that the cell volume was in accordance with the literature. From the 3D model, a cross sectional area was taken, in order to obtain the 2D model. This model was used in both first and third geometries. For the second geometry, a "parachute" shape was adopted, due to orientation of the erythrocyte. This choice was made due to the fact that erythrocytes are known to deform under flow (LIU; LIU, 2006). Moreover, it has been stated that, especially in microvessels, erythrocytes tend to form a "parachute" shape (NOGUCHI; GOMPPER, 2005). All the geometrical data is available in Table 1: In order to solve conservation equations, numerical meshes for each scenario were developed. In this step, the spatial domain is divided into infinitesimal elements, in which conservation equations are solved for each element. For each case in Table 1, three meshes were created: a coarse, an intermediate and a refined mesh. The choice of creating three meshes was to ensure that the GCI (Grid Convergence Index) method could be applied (CELIK, 2008), assuring that the influence of the numerical meshes over the numerical results is negligible. The number of elements used in each mesh, as well as mesh dimensions are presented in Tables 2 and 3:  Additionally, Figures 1-3 show the numerical meshes developed:  The number of elements chosen for each mesh was based in the factor of 1.3 proposed by Celik (2008). The parameter analyzed for the GCI development was the velocity profile in the Y direction. In all three cases, errors did not surpass 2%, which indicates that the numerical grid choice will have a considerable low influence on numerical results. Thus, for all three cases the refined grid was chosen, since the computational time to run simulations was similar for all meshes. All steps of the pre-processing were carried using the Fluent v.19 software.
For the setup used in the simulations, a transient single-phase flow, in which mass and momentum conservation equations need to be solved, respectively (Eq.1-2): (∇ · u) = 0 (Eq.1) Furthermore, since plasma is the fluid of interest, a Newtonian rheological model could be adopted. Thus, the stress tensor from Eq. 2 can be written as (Eq. 3):

τ=-μ [∇ u + (∇ u) T ] -2/3 (∇ • u) I (Eq.3)
Since blood flows have a tendency of presenting low Reynolds numbers (LIU; LIU, 2006), a laminar approach was adopted. For the processing step, plasma viscosity and density of 0.0011 Pa.s and 1020 kg/m 3 were adopted (OU et al., 2016) and an inlet velocity of 0.5 m/s was used, to represent the average blood velocity. As side boundary conditions, a zero-shear boundary was adopted, since this control volume is taken to be located at a point inside the vessel, which is far from the vessel walls. Although the blood cells were taken as being solid, the same zero-shear condition was adopted on the walls of the erythrocyte, because the cells should be moving with the plasma. Thus, the velocity around them should not equal to zero. Average blood pressure of 100 mmHg was used as exit boundary condition. The simulations used a timestep of 10 -6 , due to the fact that the Courant number was highly dependent on the timestep size. This timestep ensured Courant numbers of approximately 1. All simulations were carried for a total flow time of 0.2 s, since convergence was achieved after 0.02 s of simulation. Moreover, the residence time inside the control volume is relatively small, according to Table 2, justifying the simulation time choice. Moreover, many simplifications were adopted for the simulation in this study. Inlet boundary condition was taken as being a constant velocity, whereas regular blood flow behaves in a pulsatile manner. The erythrocytes in the simulations were taken as being stationary, while in real scenarios red blood cells have a relative velocity to the plasma. Lastly, spatial domain was taken to be two dimensional, in order to facilitate the simulations.

RESULTS AND DISCUSSION
Results obtained in the three simulations are presented in Table 3 and Figures 4 to 7. From Table 3, it can be seen that the average Reynolds number for the spatial domain is relatively low, which classifies the flow as being laminar. Although Reynolds numbers of this range are unusual for fluids with densities close to water density, these results are in accordance with experimental data from literature. Noguchi and Gomper (2005) stated that experimental Reynolds numbers for erythrocyte clusters are of 10 -2 range in microvessels. Given that the volume control is in the micron range, its flow can be compared to flow in capillaries. Additionally, since the Reynolds number equation (Eq. 4) has a characteristic length in the order of micrometers multiplying the numerator, low Reynolds numbers were expected from the simulations.
Where ρf is the fluid density, μf is the fluid viscosity, uf is the fluid velocity and Lc is the characteristic length. For this case, the fluid of interest is the plasma.    The velocity contours are shown in Figures 4 to 6. It is possible to see that the erythrocyte geometry has a significant effect on the velocity magnitude. From the velocity profiles in Figure 7, a relation of geometry to velocity can be drawn. It is visible that depending on the geometry, the deacceleration of the plasma is more pronounced. In cases two and three, the geometry promotes a higher deacceleration of the plasma at the middle line. Additionally, it promotes a higher velocity around the erythrocytes in the free flow regions, away from the center line.   Streamlines are presented in Figures 8 to 10. Although cases one and two are similar, case three showed a recirculation after the cells, despite having considerably low Reynolds numbers. This recirculation pattern can be attributed to the geometry chosen for this simulation. Since the spacing between each erythrocyte is relatively narrow, ranging from one to three erythrocyte diameters, plasma might experience an expansion after leaving the erythrocyte cluster. This result could also be associated with the fact that the cells being simulated are rigid bodies, which is not the mechanical behavior of erythrocytes.

FINAL REMARKS
The present study aimed to simulate and analyze the plasma flow around erythrocytes using the computational fluid dynamic technique. From the simulations carried out, it was possible to see that the Reynolds number in the spatial domain adopted have relative low values, indicating that the flow is laminar. Such results were in accordance with literature. Furthermore, from the velocity contours and streamlines, it was noted that erythrocyte geometry had affects in the flow field. It was observed that the fluid deaccelerated in the middle line, as expected, and accelerated on the free flow regions. Although the present study is a preliminary study, this work was used to understand the physics of this type of flow and further studies are being carried.