Design efforts in blood pumps mainly focus on hydraulic efficiency, minimization of blood trauma, and thrombus formation. It is well known that the latter effects depend upon mechanical stress loading in high shear regions of the pump and secondary flow regions, respectively. Many research groups have employed and presented flow visualization techniques to study and improve their devices, focusing mainly on centrifugal pumps ^{1–17} and partly on axial pumps. ^{18–23} However, experimental flow investigations in rotary pumps, especially inside the vanes, are considerably limited because of high rotational speed and the size of these devices, providing one reason for the application of computational fluid dynamics (CFD). In addition, information about shear loading and regions of potential thrombus formation can be obtained by computational flow predictions.

Numerical methods for flow investigations concerning rotary blood pumps have been presented by several authors who used laminar calculations ^{5, 12, 14, 15, 18, 26} and standard k-ε turbulence models, ^{4, 9, 14, 24, 25} respectively. ^{27–29} The results of these investigations were compared with experimental data in a few cases ^{14, 16, 26, 28} and showed only a rough correspondence between numerical prediction and experiment, if any at all. Results of flow visualization analysis with blood pumps have been presented by several authors. ^{7, 8, 11, 16, 17, 19–21, 23, 26–28} Laser and high speed camera based 2-D image acquisition techniques have been found suitable for inflow and outflow pattern analysis. ^{20, 23} Results of these techniques, used as a tool for validating CFD computations, are shown elsewhere. ^{16, 17}

The aim of this study is to create a refined and iteratively improved computational model of a microaxial blood pump, capable of an accurate prediction of the flow situation within an experimental pump mock loop. Flow visualization techniques (PTV) are applied for verification of the computational model. Additional comments are made on the usefulness of a standard turbulent model with logarithmic wall function regarding micropumps with relatively low Reynolds numbers (Re ≈ 10^{4}).

## Materials and Methods

### Pump Model and Mock Loop

The investigated device is an axial blood pump for left ventricular application (Intracardiac Pump/LV, Impella AG, Aachen, Germany) developed at the Helmholtz Institute, ^{29–32} as shown in Figure 1. The outer diameter of the pump housing is 6.4 mm, the inflow cannula length is 60 mm, and the integrated motor drives the rotor at a maximum speed of 30,000 rpm.

For the purpose of this study, a 1:1 scale shaft-driven replica of the original impeller was integrated into a hydraulic mock loop, shown in Figure 2. The main parts of the loop are the test compartment (Figure 3) connected to a fluid reservoir, a servomotor driven throttle valve to control afterload, an ultrasonic flowmeter (T110, Transonic Systems Inc., Ithaca, NY), pressure sensors (PVB GmbH, Kirchseeon, Germany), and a personal computer based control unit. The original micromotor is replaced by an external driving unit that consists of an electronically commutated brushless direct current (DC) motor (Faulhaber GmbH, Schönaich, Germany) located outside the test loop. Torque is transferred to the inner drive shaft by magnetic coupling. This specific design was chosen to avoid any sealed bearings on the drive shaft and to provide easy access to the motor unit. The inner shaft is supported by two bearings and covered by a cylinder made of stainless steel (AISI 301), forming a mock-up of the original motor housing. The original inflow cannula and rotor housing are replaced by a translucent block made of polymethyl methacrylate (PMMA), with an axially placed inflow boring. The outflow diameter is 25.4 mm, which approximately simulates the aortic diameter. Pressure bores are located at two inflow and two outflow positions. The fluid reservoir contains the water-glycerol test fluid, with a density of 1,059 kg/m^{3} and a viscosity of 3.6 cps, at a temperature of 21°C. To maintain a constant viscosity during pump operation, evaporation from the reservoir is prevented with a sealed lid and the fluid temperature is controlled by cooler elements within the large reservoir. Temperature is continuously measured by a standard thermoelectric thermometer. All devices were connected to a standard PC data logging unit via an analog-to-digital (AD/DA) converter (AT-MIO 16 F5, National Instruments, Austin, Texas). A controller console based on the Labview data acquisition platform (Labview 5.1, National Instruments) was established to automatically monitor pressure sensors, flow probe, and voltage controlled oscillator (VCO) motor signal and to control motor speed and valve position. In automatic operation mode, complete HQ characteristics can be measured without user interaction. HQ characteristics were taken at a rotational speed of 30,000 rpm, which is the standard operating speed of the pump *in vivo* and at 25,000 rpm and 32,500 rpm, additionally. Data logging was started with the throttle valve in the fully open position. The valve was then continuously shut down, and the resulting pressure head and flow rate were recorded.

### Visualization

Particle tracking velocimetry (PTV) is used for quantitative investigation of flow patterns and exit angles in the outlet region, and to detect the critical swirl patterns in the inflow cannula. Particle tracking velocimetry was performed using Sephacryl S-1,000 tracer particles (Amersham Biotech, Little Chalfont, Buckinghamshire, United Kingdom) 20 to 50 microns in size, with a density of 1,000 kg/m^{3}. 3-D image acquisition was established with a combination of a charge coupled device (CCD)-shutter-camera (SSC-M370, Sony, Tokyo, Tokyo, Japan) (alternatively video-shutter-camera [EX1, Canon, Tokyo, Japan]) and a mirror set-up, to vertically project the side view of the test compartment. Images were taken at a frame rate of 25/sec with a shutter aperture time of 1 to 4 msec, depending upon the flow rate, providing a track length of 3 to 10 mm. The images were recorded using a standard video cassette recorder (VCR) with freeze frame capability (FS200, Panasonic, Sagan, Japan) and then transferred to a PC with a frame grabbing extension (PC-EYE2, Eltec, Mainz, Germany). The relative motion between blades and particles could not be resolved with this method because of the high rotational speed. The images of corresponding vertical and horizontal planes of the test compartment were then spatially calibrated, and the velocity was calculated from the resulting track length in 3-D coordinates. The mirror set-up allowed the use of only one camera to avoid synchronization problems. Because of the relatively small difference between the refraction index of PMMA and water-glycerol, the distortion of the acquired camera image can be neglected.

### Computational Fluid Dynamics

For grid generation, boundary set-up, solving, and postprocessing the TASCflow package (AEA, Otterfing, Germany) is used. Pre– and postprocessing are performed on OGL workstations (Indigo2, Silicon Graphics Inc., Silver Spring, Alabama, USA) and solving is carried out on clustered workstations (>1,000 MFlops) at the computing center of Aachen University of Technology. The computational domain, shown in Figure 4, includes the vane regions, the inflow cannula, and the outlet region including the outlet of the test compartment. The inclusion of the test compartment into the computational domain entails a considerable computational expenditure compared with the actual pump itself. There are two reasons for extending the domain to the boundaries of the entire test compartment. First, to ensure proper validation, the computational domain has to match the real geometry as closely as possible. This includes an exact match between actual locations of pressure sensors, which are of course upstream and downstream of the actual pump, and their computational counterparts. Second, an accurate outflow boundary condition requires a region where the assumption of a constant pressure distribution is fully valid. This is most probably the case well downstream of the blades and the outlet swirl at the outlet tube of the test compartment. The grid domain is reduced to a 180° section of the pump to take advantage of periodicity. This condition requires a size reduction of the modeled outlet tube, compared with the real compartment, to take into account the fact that actually two outlet tubes are modeled because of periodicity.

The entire mesh, as shown in Figure 4, is H-shaped, with an O-grid topology around the blade. It is composed of 10 blocks in total, incorporating a hexahedral structure with refinements toward near wall regions, especially in blade regions. The number of nodes is iteratively increased to a maximum of 310,000 until the numeric solution is found independent of the grid density. Major effort was put into designing a grid with skew angles not less than 30°, with a majority of 80% being in the range of 65° to 90° as convergence and stability are found to depend upon the grid quality in this respect. In addition, because of implementation of the wall functions, the angle between element surfaces along and adjacent to the walls were kept close to 90°.

The 0.1 mm gap between blade tip and housing is modeled by 9 nodes. Boundary conditions are defined as follows: velocity distribution is set constant at the entrance of the inflow cannula. Eddy length of 10^{−5} m and turbulent intensity of 10^{−3} % is set at the inflow boundary region. Outlet pressure is set constant in the outlet tube of the test compartment. The vane section of the pump is defined as a constantly rotating frame to facilitate steady flow computation. In the subsequent solving process the continuous governing differential equations are replaced by their discrete counterparts and discretized with the Finite Volume method and a Finite Element approach for spatial discretization. The corresponding stationary governing equations in vector notification for incompressible flow (ρ = const.) with constant viscosity (η = const.) are given below. MATH

In the conservation of momentum equation (Equation 2) the last term denotes the influence of the Reynolds stress tensor ς_{∼turb}, which is the result of the time averaged formulation of turbulent fluctuation. It only applies in turbulent calculations, however. For laminar flow this term becomes zero and equation 2 forms the incompressible Navier-Stokes equation. In this study, a k-ε model is used to relate the Reynolds stress terms with the mean flow variables by using an eddy viscosity assumption. However, some comments have to be made on laminar *versus* turbulent calculations. Generally, a turbulent calculation is only appropriate in flow conditions with a distinct (generally turbulent) boundary layer and fully developed turbulent flow regime away from the wall. To fully resolve the solution in the near-wall region is very costly, because the required number of nodes would be quite large. Thus, a common approach is to model this region using wall functions, where the near wall tangential velocity is related to the wall shear stress by means of a logarithmic relationship, or by employing a modified turbulence model in near wall regions that can more closely resolve the flow within the boundary layer (two-layer model). Laminar calculation, however, does not involve any resolution of turbulent fluctuations and is therefore only appropriate in flow conditions where the turbulence is adequately low. On the other hand, laminar calculations used in conjunction with fine grids have a much better capability of resolving the boundary layer if viscous forces are dominating the turbulent fluctuations. In former studies some blood pump simulations have been performed based on a laminar flow approach, ^{12, 14, 18} and some on a turbulent (k-ε) model. ^{9, 24} Now, which approach is suitable for the calculation of the discussed micropump?

The machine Reynolds number, based on the running length of the blades and the average meridional speed, which can be taken as a rough characteristic for the flow regime, is fairly low—between 5,000 and 15,000, depending on the operating point, which could justify a laminar calculation. This is because most turbulence models are more suitable for high Reynolds flows and thus laminar calculation might offer some advantage even if the flow is not laminar. Another important criterion, the critical running length where laminar/turbulent transition takes place, can be roughly estimated based on equations for the flat plate:MATH where *u* denotes the tangential velocity far away from the plate, *x*, the running length, and *v*, the dynamic viscosity. In this case, the critical length results to tenths of meters (approximately 0.2 m), or assuming the flat plate is an appropriate model for a pump blade, the boundary layer at the blades is not turbulent, but laminar. On the other hand, the boundary layer is very thin compared with the vane dimensions. Approximations again based on the flat plate yield:MATHwhere δ_{99} denotes the thickness of the (laminar) boundary layer at a point where the tangential velocity reaches 99% of the velocity, far away from the wall. Thus, the grid resolution of the wall regions would have to be extraordinarily high. Furthermore, the circumstance of the entire test compartment being modeled has another drawback; the discontinuous diffuser-type opening, downstream of the pump, entails a large increase in turbulence and high vorticity in combination with a flow separation. This flow region is not suitable to be resolved with laminar calculation. Within this study, convergent solutions on the entire test compartment could only be accomplished with the use of turbulence models.

A two layer model, with a one equation wall model, ^{30} is used for turbulent calculations in this study, and the standard k-ε model is used in regions away from the wall. With this approach, the influence of wall effects is taken into consideration, as well as the highly turbulent flow field in the pump outlet. Discretization of the previous governing equations was performed using a mass weighted differencing scheme (MWS), blended with a 5% upwind differencing scheme (UDS), for stability reasons. Both schemes operate on a finite volume basis.

## Results

At a given rotational speed, the pressure head is a function of the flow rate. This function is represented in pump specific HQ characteristic curves, which have been measured and computed in this study for rotational speeds of 25,000 rpm, 30,000 rpm, and 32,500 rpm. Figure 5 shows a dimensionless representation of the results. Measured curves are represented by a continuous line, and CFD calculated operating points are marked as dots. The dimensionless flow coefficient (abscissa) was calculated from the following equation:

where *c* denotes the relative meridional flow velocity, which can be obtained from MATH

*u*_{B} represents the rotational speed of the blade at a diameter of *D*_{B}, and *Q* represents the flow rate. The pressure head coefficient, which is shown on the ordinate of Figure 5, was calculated from the following equation:

with Δp as the pressure head between inlet and outlet measurement locations.

The plot of measured and computed HQ data demonstrates good agreement; the characteristics of the HQ curve are perfectly represented by CFD calculations. However, there is a notable difference in magnitude of the flow coefficient that presumably has two main reasons. (1) The two layer representation of the turbulent flow in the vane regions implies certain limitations in accuracy at low Reynolds numbers. Viscous effects on the walls are underestimated, which results in lower total pressure gain. (2) The outlet region of the impeller is highly turbulent (Reynolds number of approximately 50,000–100,000), and the discontinuous diffuser flow region in the vicinity of the trailing edge leads to vorticity in the outlet. This region presumably causes some flaws in the CFD prediction.

To further compare the validity of the computational model, PTV images, taken from the outlet region of the pump, are compared with streak lines generated from the CFD data in Figure 6. On the upper side, a PTV image of a particle exiting the pump vane is shown. In comparison, the lower image shows a sample of streak lines computed from the absolute speed vectors within the vane region based on CFD. The PTV images represent a single particle track, disregarding radial coordinates, and can only be taken as a rough characterization of the outlet flow direction because of inertia of the particles. The images in Figure 6 represent a median particle track, from 30 PTV images that were taken per operation point. Operation point for both images is 4.5 L/min.

The comparison of computed exit angles and measured pathlines in Figure 6 demonstrates an acceptable agreement. The outlet angles determined by the PTV method have shown an estimated deviation of approximately 5°, whereas CFD (past validation) allows the specific determination of angles at any given radius. Exemplifying operating point of 4.5 L/min as shown in Figure 6, an angle of 42° ± 3° was determined using 30 PTV images. Detailed CFD results at the same operating point are shown in Figure 7. Angles at 20 meridional span coordinates were calculated from CFD data (continuous line with markers). As a reference, the upper and lower limits of the range of values that were determined by PTV are also shown in Figure 7 (horizontal lines). The CFD results show that the outflow region near the hub is affected by wall friction. The peak at a span value of 0.1 in Figure 7, represents the resulting increase in circumferential velocity. In the remaining vane region the angle is gradually increasing, with an increasing radial coordinate, which demonstrates the influence of circumferential velocity ω · r on the outflow angle.

The presented micropump was designed to ensure inflow conditions without an angular momentum (swirl) when the flow hits the leading edge of the blades. However, this condition changes at certain pressure levels. It is important to verify at which operating points the inflow swirl starts. Results of this study show that the used CFD model is capable of predicting these critical operating points quite accurately. Figure 8 demonstrates the appearance of inlet swirl at an operating point of 1 L/min. A PTV image and a CFD streak line plot illustrate the change in angular momentum at approximately 1 diameter upstream of the leading edge of the blades.

Figure 9 shows the pressure rise through the vanes of the micropump. It can be seen that the minimum pressure occurs at the leading edge of the blades because of acceleration and is then increasing through the vanes, and ends with a maximum gainable pressure rise in the vane at the trailing edge of the pump. These results are very important for blade design considerations and tracking energy losses of the entire system.

The overall efficiency, which relates the actual specific energy rise of the entire turbo machine to the mechanical power on the shaft, can be written as:MATH

where the following notations are used:

- Y
_{t}: specific energy rise, including pressure head and absolute speed
- ṁ mass flow
- P
_{mech}: mechanical power on shaft
- Δp: pressure head (non-specific)
- Δc: change of absolute speed between inlet and outlet of the pump
- ω: angular velocity
- ρ: density
- M: torque on shaft

Much more commonly used is the hydraulic efficiency, which relates the hydraulic energy to the shaft power according to the following equation:MATH

In connection with findings of hydraulic efficiency being a necessity, but not a sufficient and sensitive criterion for hemolysis minimization, ^{13, 34} an aim of pump design needs to be a configuration with maximum efficiency at the standard operating point of the pump. In addition, the efficiency should not vary too much when the operation point is changed. Although determination of the energy rise can be measured, the instrumental determination of mechanical shaft power or torque is extremely costly and inaccurate with pumps operating at several thousand rpm. In this study, CFD data are used to calculate efficiency from power values for different operating points.

For the discussed micropump, this efficiency is shown in Figure 10, lower graph, as a function of the flow in nondimensional notation. The graph shows a weak maximum at a flow of approximately 0.18, which corresponds to 3.5 liters/min. The graph has a reasonably flat character, which indicates a more or less constant efficiency for a range of operating points. This is a good precondition for *in vivo* pump operation. The upper graph in Figure 10 shows the ratio between hydraulic energy and the Euler specific energy, which is derived from the following equation:MATH

where H is the Euler specific energy rise; u is the rotational speed (ωr); and c_{u} is the circumferential component of the absolute speed vector in inlet (Index 1) and outlet (index 2).

Figure 10 shows the latter efficiency again as a function of the nondimensional flow. The graph shows a constant efficiency of approximately 80%, regardless of flow rate.

## Conclusions

The aim of this study was the set-up and experimental verification of a computational model capable of accurate flow prediction in a microaxial blood pump. The main focus was set on the validation of the developed model, because CFD predictions are not reliable without additional confirmation. It has been found that it is indispensable for proper validation to perfectly match the computational domain, including the up- and downstream periphery of the test environment. It has also been shown that flow visualization, in this case PTV, is generally suitable for detailed verification of even small devices like the micropump discussed.

The necessity of careful consideration regarding turbulence modeling and wall function solutions has also been discussed. The results demonstrate very good agreement between calculation and measurements. With the presented set-up, it has not only become possible to predict HQ characteristics quite accurately, but also to get detailed information about flow patterns and important numbers concerning efficiency. However, the total effort to ensure a reasonably accurate flow prediction by means of CFD is quite high. This set-up can and will now be used for detailed analysis involving flow optimization and basic approaches for the prediction of shear induced hemolysis.

## References

1. Takami Y, Yamane S, Makinouchi K, Glueck J, Nose Y: Mechanical white blood cell damage in rotary blood pumps. Artif Organs 21: 138–142, 1997.

2. Andrade A, Biscegli J, Sousa JE, Ohashi Y, Nose Y: Flow visualization studies to improve the spiral pump design. Artif Organs 21: 680–685, 1997.

3. Nakamura S, Ding W, Smith WA, Golding LA: Numeric flow simulation for an innovative ventricular assist system secondary impeller. ASAIO J 45: 74–78, 1999.

4. Allaire PE, Wood HG, Awad RS, Olsen DB: Blood flow in a continuous flow ventricular assist device. Artif Organs 23: 769–773, 1999.

5. Wood HG, Anderson J, Allaire PE, McDaniel JC, Bearnson G: Numerical solution for blood flow in a centrifugal ventricular assist device. Int J Artif Organs 22: 827–836, 1999.

6. Ahmed S, Funakubo A, Sakuma I, Fukui Y, Dohi T: Experimental study on hemolysis in centrifugal blood pumps: improvement of flow visualization method. Artif Organs 23: 542–546, 1999.

7. Araki K, Taenaka Y, Masuzawa T, et al: A flow visualization study of centrifugal blood pumps developed for long-term usage. Artif Organs 17: 307–312, 1993.

8. Ikeda T, Yamane T, Orita T, Tateishi T: A quantitative visualization study of flow in a scaled-up model of a centrifugal blood pump. Artif Organs 20: 132–138, 1996.

9. Miyazoe Y, Sawairi T, Ito K, et al: Computational fluid dynamic analyses to establish design process of centrifugal blood pumps. Artif Organs 22: 381–385, 1998.

10. Mizuguchi K, Damm G, Benkowsky R, et al: Development of an axial flow ventricular assist device: in vitro and in vivo evaluation. Artif Organs 19: 653–659, 1995.

11. Mulder MM, Hansen AC, Mohammad SF, Olsen DB: In vitro investigation of the St. Jude Medical Isoflow centrifugal pump: Flow visualization and hemolysis studies. Artif Organs 21: 947–953, 1997.

12. Pinotti M, Rosa ES: Computational prediction of hemolysis in a centrifugal ventricular assist device. Artif Organs 19: 267–273, 1995.

13. Schima H, Muller MR, Papantonis D, et al: Minimization of hemolysis in centrifugal blood pumps: influence of different geometries. Int J Artif Organs 16: 521–529, 1993.

14. Takiura K, Masuzawa T, Endo S, et al: Development of design methods of a centrifugal blood pump with in vitro tests, flow visualization, and computational fluid dynamics: results in hemolysis tests. Artif Organs 22: 393–398, 1998.

15. Thomas DC, Butler KC, Taylor LP, et al: Continued development of the Nimbus/University of Pittsburgh (UOP) axial flow left ventricular assist system. ASAIO J 43: M564–M566, 1997.

16. Yamane T, Asztalos B, Nishida M, et al: Flow visualization as a complementary tool to hemolysis testing in the development of centrifugal blood pumps. Artif Organs 22: 375–380, 1998.

17. Treichler J, Rosenow SE, Damm G, et al: A fluid dynamic analysis of a rotary blood pump for design improvement. Artif Organs 17: 797–808, 1993.

18. Antaki JF, Ghattas O, Burgreen GW, He B: Computational flow optimization of rotary blood pump components. Artif Organs 19: 608–615, 1995.

19. Kerrigan JP, Shaffer FD, Maher TR, Dennis TJ, Borovetz HS, Antaki JF: Fluorescent image tracking velocimetry of the Nimbus AxiPump. ASAIO J 39: M639–M643, 1993.

20. Kerrigan JP, Yamazaki K, Meyer RK, et al: High-resolution fluorescent particle-tracking flow visualization within an intraventricular axial flow left ventricular assist device. Artif Organs 20: 534–540, 1996.

21. Butler K, Thomas D, Antaki J, et al: Development of the Nimbus/Pittsburgh axial flow left ventricular assist system. Artif Organs 21: 602–610, 1997.

22. Butler KC, Maher TR, Borovetz HS, et al: Development of an axial flow blood pump LVAS. ASAIO J 38: M296–M300, 1992.

23. Wernicke JT, Meier D, Mizuguchi K, et al: A fluid dynamic analysis using flow visualization of the Baylor/NASA implantable axial flow blood pump for design improvement. Artif Organs 19: 161–177, 1995.

24. Bludszuweit C: Model for a general mechanical blood damage prediction. Artif Organs 19: 583–589, 1995.

25. Bludszuweit C: Three-dimensional numerical prediction of stress loading of blood particles in a centrifugal pump. Artif Organs 19: 590–596, 1995.

26. Konig CS, Clark C, Mokhtarzadeh-Dehghan MR: Investigation of unsteady flow in a model of a ventricular assist device by numerical modelling and comparison with experiment. Med Eng Phys 21: 53–64, 1999.

27. Woodard JC, Shaffer FD, Schaub RD, Lund LW, Borovetz HS: Optimal management of a ventricular assist system. Contribution of flow visualization studies. ASAIO J 38: M216–M219, 1992.

28. Sturm C, Li W, Woodard JC, Hwang NH: Fluid mechanics of left ventricular assist system outflow housings. ASAIO J 38: M225–M227, 1992.

29. Meyns B, Siess T, Nishimura Y, et al: Miniaturized implantable rotary blood pump in atrial-aortic position supports and unloads the failing heart. Cardiovasc Surg 6: 288–295, 1998.

30. Patel VC, Rodi W, Scheuerer G: Turbulence models for near-wall and low Reynolds number flows: A review. AIAA J 23: 1308–1319, 1985.

31. Siess T, Reul H, Rau G: Concept, realization, and first in vitro testing of an intraarterial microaxial blood pump. Artif Organs 19: 644–652, 1995.

32. Siess T, Reul H, Rau G: Hydraulic refinement of an intraarterial microaxial blood pump. Int J Artif Organs 18: 273–285, 1995.

33. Schlichting H: Grenzschichttheorie. Springer, Berlin, 1997.

34. Bluestein M, Mockros LF: Hemolytic effects of energy dissipation in flowing blood. Med Biol Eng 7: 1–16, 1969.