Received: 8 October 2022 Accepted: 10 January 2023 DOI: 10.1002/pamm.202200295 Analysing the bone cement flow in the injection apparatus during vertebroplasty Zubin Trivedi1,∗, Dominic Gehweiler2, Jacek K. Wychowaniec2, Tim Ricken3, Boyko Gueorguiev-Rüegg2, Arndt Wagner4,5, and Oliver Röhrle1,5 1 Institute for Modelling and Simulation of Biomechanical Systems, University of Stuttgart, 70569 Stuttgart 2 AO Research Institute (ARI), Clavadelerstrasse 8, 7270 Davos, Switzerland 3 Institute of Structural Mechanics and Dynamics of Aerospace Structures, University of Stuttgart, Pfaffenwaldring 27, 70569 Stuttgart 4 Institute of Applied Mechanics (CE), University of Stuttgart, Pfaffenwaldring 7, 70569 Stuttgart 5 Stuttgart Center for Simulation Science (SC SimTech), Pfaffenwaldring 5a, 70569 Stuttgart Vertebroplasty, a medical procedure for treating vertebral fractures, requires medical practitioners to inject bone cement inside the vertebra using a cannula attached to a syringe. The required injection force must be small enough for the practitioner to apply it by hand while remaining stable for a controlled injection. Several factors could make the injection force unintuitive for the practitioners, one of them being the non-Newtonian nature of the bone cement. The viscosity of the bone cement varies as it flows through the different parts of the injection apparatus and the porous cancellous interior of the vertebra. Therefore, it is important to study the flow of bone cement through these parts. This work is a preliminary study on the flow of bone cement through the injection apparatus. Firstly, we obtained the rheological parameters for the power law model of bone cement using experiments using standard clinical equipment. These parameters were then used to obtain the shear rate, viscosity, and velocity profiles of the bone cement flow through the cannula. Lastly, an analysis was carried out to understand the influence of various geometrical parameters of the injection apparatus, in which the radius of the cannula was found to be the most influential parameter. © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. 1 Introduction Vertebral fractures are a common occurrence worldwide, generally caused by accidents or diseases like osteoporosis. Verte- broplasty is a medical procedure for treatment and prevention of vertebral fractures. Vertebroplasty entails injecting a “bone cement" into the porous cancellous part of the vertebra. This bone cement is a curing polymer and is left to set inside the vertebra. Once it is solid, it restores the structural integrity of the vertebra. Details about the procedure can be found in [1]. Improper injection of the bone cement could lead to leakage into unwanted places, causing complications like embolism and paralysis [2, 3]. Of the many factors that make the procedure challenging for practitioners, one of the most important is the non-Newtonian rheology of bone cement. This means that the viscosity has a non-linear relationship with the shear rate. The shear rate depends on the flow rate and the geometry of the structure through which the bone cement flows. The resultant viscosity affects the required injection force and the haptic feedback to the practitioner. Therefore, it is important to understand the flow of bone cement in all the parts that it flows through. The bone cement flow path can be divided into two parts: (i) through the porous cancellous bone; (ii) through the injection apparatus, i.e. the syringe, and the cannula. The cancellous part of the vertebra is a network of trabeculae with very small pores filled with bone marrow. Naturally, the flow of the bone cement in such a structure is dependent upon the non-Newtonian characteristics, the geometry of the structure, and the consistency of the bone marrow. Such a complex problem could be addressed using a continuum porous media flow model, like in [4] using the Theory of Porous Media [5]. This is a macro-scale model, hence the Darcy equations are used to govern the bone cement flow. Due to the high viscosity and the injection pressure, the capillary pressure can be neglected, giving an advection-dominated flow. The numerical treatment of such a flow needs to be done with specific regard to the choice of variables and the discretization technique. Furthermore, the viscosity at the pore-scale needs to be upscaled to the macro-scale, which is not trivial for the case of non-Newtonian fluid [6]. On the other hand, the flow of the bone cement through the injection apparatus is not extensively studied, although it could account for up to 85% of the required injection pressure [7], thus making it a significant aspect to determine the injection force required by the vertebroplasty practitioner. Such a study would aid the decision-making in the selection and design of injection apparatus [8]. This work aims to analyse and understand the bone cement flow through the injection apparatus, specifically the syringe, and the cannula; using a power law model. For this purpose, the rheological parameters are obtained by pressure measurements from experiments at various flow rates. The simple cylindrical geometry of the syringe and the cannula allows the analytical calculation for the required injection pressure and the profiles of the shear rate, viscosity, and the velocity along the radius of the cylinder. The analytical calculations are then used to study the profiles of rheological parameters and also to determine the most influential parameter in the geometry of the injection apparatus. ∗ Corresponding author: e-mail zubin.trivedi@imsb.uni-stuttgart.de, phone +49 711 685 60948, fax +49 711 685 66347 This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. PAMM · Proc. Appl. Math. Mech. 2022;23:1 e202200295. www.gamm-proceedings.com 1 of 6 https://doi.org/10.1002/pamm.202200295 © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. 2 of 6 Section 2: Biomechanics 2 Methods 2.1 Experiments (a) (b) Fig. 1: (a) The injection apparatus used in the experiment and its simplified form used for analytical calculations (b) The custom-made injector setup used for injecting the bone cement at constant flow rate For determining the rheological parameters of the bone cement, we measured the injection force required for injecting the bone cement through the injection apparatus as shown in Figure 1a at two different flow rates. The cannula used here was an 8-gauge cannula. For this, we used the experimental setup shown in Figure 1b. In this setup, we designed and manufactured a custom-made bone cement injector. The injector consisted of a carriage driven by a stepper motor using a ball screw with a 5 mm pitch. The stepper motor gave 200 steps per revolution, such that the carriage could be moved with a precision of 0.025 mm. A 200 N load cell measured the force applied to the syringe. The pressure was then calculated by dividing the measured force by the cross-sectional area of the syringe. The bone cement Vertecem V+ was used for the experiment. The cement was prepared by mixing 2.6 g of the polymer powder with 1 ml of the monomer liquid, which were provided in the cement preparation kit by DePuy Synthes. This mixing method is different from the manufacturer-prescribed way using the kit, but we used our mixing method to save material and verified by rheometer measurements that it consistently reproduces cement behaviour as is intended in clinical application. The time from the start of mixing to the start of injection was measured using a stopwatch, which was about 200-210 seconds. The procedure was repeated to carry out two experiments with injections at different flow rates: one at a flow rate of 0.1 ml/s and the other at a flow rate of 0.4 ml/s. 2.2 Analytical calculations For the analytical calculation of the bone cement flow, the geometry of the apparatus was simplified to the system as shown in Figure 1a. The syringe and the cannula are assumed to be cylinders with the given dimensions. The cylindrical geometry of the cannula in the simplified setup includes the cannula in the original setup, and the coupling between the syringe and the cannula, since the flow path inside the coupling is also a cylinder with same radius as that of the cannula. The nozzle of the syringe is short in length and hence ignored. The change of cross-section from syringe to cannula would add to the pressure loss in reality, however this has not been included in this study. Although the cannula had a side-opening, we assume that the cannula is open-ended. Fig. 2: The diagram shows the flow of a fluid through a cylindrical channel of radius R. The analytical calculation is done over a control volume of length L. © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. www.gamm-proceedings.com PAMM · Proc. Appl. Math. Mech. 23:1 (2022) 3 of 6 As shown in Figure 2, we assume that the fluid flow throws through a cylinder with radius R. A control volume of length L is taken where the left end is subjected to higher pressure p+∆p compared to the right end with pressure p, resulting in a flow from left to right. The pressure is assumed to change only along the length of the cylinder. The fluid is assumed to have no slip at the walls. This gives a velocity profile along the radius with maximum velocity at the centre and zero at the walls. The viscosity of the fluid relates the viscous shear stress to the shear rate. This is written as τ = µγ̇. (1) where τ is the viscous shear stress, µ is the viscosity, and γ̇ is the shear rate. Since the bone cement is a non-Newtonian fluid, the viscosity is not constant, rather it is dependent on the shear rate. The dependence of viscosity on the shear rate is described by a rheological model. For the selection of the rheological model for bone cement, the power law model is an appropriate choice since it has already been validated for bone cement [10]. The equation of the power law model is given as µ = Kγ̇n−1. (2) where K is the viscosity at unit shear rate, and n is the power law index. Substituting the power law model in Equation 3, we get τ = Kγ̇n. (3) Now, for any cylinder of radius r < R concentric to the control volume, the resultant pressure acting on the cross-sectional area needs to be equilibrated by the viscous shear stress caused by the velocity profile. This gives (p+∆p− p)πr2 = τ(2πrL), (4) ⇒ τ = r∆p 2L . (5) On the other hand, the shear rate can be simply obtained from the tangent of the velocity profile, using γ̇ = −dv(r) dr . (6) Hence, from Equations 3, 5 and 6, we get r∆p 2L = K ( − dv(r) dr )n , ⇒ −dv(r) dr = ( ∆p 2KL ) 1 n r 1 n . (7) Integrating Equation 7 with respect to radius r, we get v(r) = − ( ∆p 2KL ) 1 n ( n 1 + n ) r 1+n n + C. (8) where C is the integration constant. By applying the boundary condition that velocity is zero at the walls, i.e. v = 0 at r = R, we get 0 = − ( ∆p 2KL ) 1 n ( n 1 + n ) R 1+n n + C, ⇒ C = ( ∆p 2KL ) 1 n ( n 1 + n ) R 1+n n (9) Inserting the value of C in Equation 8 and rearranging gives v(r) = ( ∆p 2KL ) 1 n ( n 1 + n ) R 1+n n ( 1− r 1+n n R 1+n n ) . (10) The velocity is maximum at the centre of the cylinder, therefore v = vmax at r = 0. Putting this in Equation 10, we get the value of the maximum velocity at the centre as vmax = ( ∆p 2KL ) 1 n ( n 1 + n ) R 1+n n . (11) www.gamm-proceedings.com © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. 4 of 6 Section 2: Biomechanics Hence, we can write the velocity at a radius r in terms of vmax, i.e. v(r) = vmax ( 1− r 1+n n R 1+n n ) . (12) The volumetric flow rate can then be obtained by integrating the product of the velocity and the area of a ring with radius r and infinitesimal thickness dr, over the entire radius of the cylinder. This is written as Q = ∫ R 0 2πrv(r)dr. (13) Substituting the velocity from Equation 12, Q = ∫ R 0 2πrvmax ( 1− r 1+n n R 1+n n ) dr, = 2πvmax ∫ R 0 ( r − r 1+2n n R 1+n n ) dr, = 2πvmax [ r2 2 − n 1 + 3n r 1+3n n R 1+n n ]R 0 , = 2πvmax [ R2 2 − n 1 + 3n R 1+3n n R 1+n n − 0 + 0 ] , = ( 1 + n 1 + 3n ) πR2vmax. (14) The average velocity vavg is obtained by dividing the flow rate by the total cross-sectional area A of the cylinder. Thus, we get vavg = Q A = Q πR2 = ( 1 + n 1 + 3n ) vmax. (15) Finally, by substituting vmax from Equation 11 in Equation 15 and rearranging, we get the pressure loss in the cylinder in terms of the flow rate, called the Rabinowitsch equation [9], given as ∆p = 2LK R ( 1 + 3n n Q πR3 )n . (16) We also wanted to look at the profiles of the shear rate, viscosity, and velocity along the radius of the cannula. Using Equations 6, 7, and 16, we obtain the shear rate profile along the radius, i.e. γ̇(r) = ( 1 + 3n n Q πR3 )( r R ) 1 n . (17) Subsequently, we get the viscosity profile by substituting the shear rate in the Equation 2, µ(r) = K ( 1 + 3n n Q πR3 )n−1( r R )n−1 n . (18) Finally, we get the velocity profile by substituting ∆p from Equation 16 into Equation 10, which upon simplification yields v(r) = 1 + 3n 1 + n Q πR2 ( 1− r 1+n n R 1+n n ) . (19) Note that these equations hold only for fully developed flow, and hence may not be true near the start of the tube. Further details about the derivations can be found in [9] and [11]. For our case, the pressure loss is the sum of the pressure losses in the syringe and the cannula. For two flow rates Q1 and Q2, we get a system of two equations, i.e. ∆p1 = ∆ps1 +∆pc1 = 2LsK Rs ( 1 + 3n n Q1 πR3 s )n + 2LcK Rc ( 1 + 3n n Q1 πR3 c )n , (20) ∆p2 = ∆ps2 +∆pc2 = 2LsK Rs ( 1 + 3n n Q2 πR3 s )n + 2LcK Rc ( 1 + 3n n Q2 πR3 c )n . (21) where Ls, Rs are the length and radius of the syringe; Lc, Rc are the length and radius of the cannula. © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. www.gamm-proceedings.com PAMM · Proc. Appl. Math. Mech. 23:1 (2022) 5 of 6 3 Results The injection pressures obtained for the different flow rates are given in Table 1. These values (after conversion to SI units) were inserted in Equations 20 and 21 along with the values of the dimensions of the syringe and the cannula. Thus, we obtained a system of two non-linear equations with two unknowns K and n. These equations were solved using an iterative solver to obtain the values K = 1554 and n = 0.30. These values could then be used for further calculations and analysis. We could also calculate the pressure loss in the syringe and the cannula separately using Equation 16. For Q = 0.1 ml/s, we obtained ∆ps = 0.12 MPa and ∆pc = 1.24 MPa. Similarly, for Q = 0.4 ml/s, we obtained ∆ps = 0.18 MPa and ∆pc = 1.86 MPa. Evidently, the pressure loss in the syringe is less than 10% of that in the cannula, meaning that most of the pressure loss is occurring in the cannula. Rs Ls Rc Lc Experiment 1 Experiment 2 Parameters (mm) (mm) (mm) (mm) Q1 (ml/s) ∆p1 (MPa) Q2 (ml/s) ∆pc2 (MPa) K n 3.2 62 1.6 195 0.1 1.4 0.4 2.1 1554 0.30 Table 1: Summary of values used in Equations 20 and 21 to obtain the rheological parameters K and n for the bone cement The shear rate, viscosity, and velocity profiles for flow inside the cannula are shown in Figure 3 for 0.1 ml/s flow rate. The shapes of these profiles do not change with flow rate, hence, they are the same for 0.4 ml/s flow rate. The shear rate profile is S-shaped, as opposed to a linear profile one observes in the case of a Newtonian fluid. The shear rate is maximum at the walls, and zero at the centre. Correspondingly, the viscosity is the least at the walls, while it tends towards infinity near the centre. This is one of the drawbacks of the power law model since there is no upper viscosity limit. Finally, the velocity profile looks like a parabola, but is not exactly parabolic. This is due to the region around the centre where the velocity is maximum and also does not vary much with the radius. a) b) c) Fig. 3: Profile along the radius for (a) Shear rate (b) Viscosity (c) Velocity for bone cement flow though the cannula at 0.1 ml/s flow rate Finally, we wanted to see how varying the geometry affects the results, and which geometric parameters have the most significant effect. Figure 4 shows the change in the injection pressure caused by varying the geometry of the syringe and the cannula setup through their lengths and radii. Each parameter is varied over two orders of magnitude. Although most of the values within this range are impractical for the given components, it helps us understand how much influence each of these parameters have. In this case, a steep curve indicates a strong influence, while a flat curve implies a weak influence. We found that the injection pressure is much more sensitive of the radii rather than the lengths, especially the radius of the cannula. The syringe has as much influence only when Rs < Rc = 1.6 mm, which is quite unrealistic. For the realistic case of Rs > Rc, the curve flattens, implying that the influence of Rs is greatly diminished. On the other hand, the length of the syringe has very little influence, as evident from its nearly flat curve. The length of the cannula has a much greater impact, which is probably because of the smaller radius of the cannula. However, the impact of the lengths is still less compared to the impact of the radii. 4 Conclusion and Outlook The analysis of the pressure loss occurring in various stages of the flow path of the bone cement is important to understand and control the injection force required by the practitioner of vertebroplasty. The injection pressure required for pushing the bone cement through the injection apparatus is not negligible and may have a significant influence on the force needed by the practitioner. The majority of the pressure loss occurs in the cannula owing to it being longer and smaller in radius. The radius of the cannula is the most influential geometric parameter affecting the injection pressure. Although this was a preliminary analysis, a more detailed study using Computational Fluid Dynamics simulations could provide more insights and help make decisions regarding the selection of size and length of the cannula and the volume of the syringe, while also including the effects of a real geometry which are difficult to include in the analytical calculation, e.g. change in cross-section between syringe and cannula, or having a side-opening in the cannula. www.gamm-proceedings.com © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. 6 of 6 Section 2: Biomechanics a) b) Fig. 4: Injection pressure for (a) varying radii of the syringe and the cannula (b) varying lengths of the syringe and the cannula. For graph corresponding to each geometric parameter, the other geometric parameters are kept constant as in the values in Table 1. Acknowledgements This work has been funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – Project Number 327154368 – SFB 1313. We also acknowledge the support by the Stuttgart Center for Simulation Science (SimTech). J.K.W. ac- knowledges the European Union’s Horizon 2020 (H2020-MSCA-IF-2019) research and innovation programme under the Marie Sklodowska- Curie grant agreement 893099 — ImmunoBioInks. Open access funding enabled and organized by Projekt DEAL. References [1] M. E. Jensen, A. J. Evans, J. M. Mathis, D. F. Kallmes, H. J. Cloft, and J. E. Dion, Am J Neuroradiol 18, 1897–1904 (1997). [2] J. Bernhard, Ann Rheum Dis 62, 85–86 (2003). [3] J. Ratliff, T. Nguyen, and J. Heiss, Spine J 26, e300–e302 (2001). [4] Z. Trivedi et al., arXiv:2209.14654, arXiv (2022). [5] W. Ehlers and J. Bluhm, Porous media: Theory, experiments and numerical applications (Springer, Berlin, 2002). [6] J. Pearson and P. Tardy, J Nonnewton Fluid Mech 102, 447–473 (2002). [7] Z. Lian, C. K. Chui, and S. H. Teoh, Comput Biol Med 38, 304–312 (2008). [8] J. Krebs, S. J. Ferguson, M. Bohner, G. Baroud, T. Steffen, and P. F. Heini, Spine 30, E118–E122 (2005). [9] R. S. Lenk, in: Polymer Rheology, (Springer, Dordrecht, 1978), 75–85. [10] G. Baroud and F. B. Yahia, Proc Inst Mech Eng H: J Eng Med 218, 331–338 (2004). [11] M. Simpson and W. Janna, Proc Int Mech Eng Congress Expo 9, 173-180 (2008). © 2023 The Authors. Proceedings in Applied Mathematics & Mechanics published by Wiley-VCH GmbH. www.gamm-proceedings.com