Background: Today, the most common method for kidney stone therapy is extracorporeal shock wave lithotripsy. Current research is a numerical simulation of kidney stone fragmentation via ultrasonic shock waves. Most numerical studies in lithotripsy have been carried out using the elasticity or energy method and neglected the dissipation phenomenon. In the current study, it is solved by not only the linear acoustics equation, but also the Westervelt acoustics equation which nonlinearity and dissipation are involved.
Objective: This study is to compare two methods for simulation of shock wave lithotripsy, clarifying the effect of shock wave profiles and stones’ material, and investigating side effects on surrounding tissues.
Material and Methods: Computational study is done using COMSOL Multiphysics, commercial software based on the finite element method. Nonlinear governing equations of acoustics, elasticity and bioheat-transfer are coupled and solved.
Results: A decrease in the rise time of shock wave leads to increase the produced acoustic pressure and enlarge focus region. The shock wave damages kidney tissues in both linear and nonlinear simulation but the damage due to high temperature is very negligible compared to the High Intensity Focused Ultrasound (HIFU).
Conclusion: Disaffiliation of wave nonlinearity causes a high incompatibility with reality. Stone’s material is an important factor, affecting the fragmentation.
Kidney stone disease is a common problem worldwide with a prevalence of approximately 10-12% of men and 5-6% of women in Western countries [ 1 ]. Several factors such as genetic, nutrition, and geographic and socioeconomic have influences on kidney stone disease [ 2 ]. There are different types of minerals that can form urinary stones and Calcium Oxalate is the most common stone component [ 3 ]. Y Warty et al. [ 4 ] analyzed the organic material and declared kidney stones have a different composition of elements in each layer.
Ultrasonic waves have a variety of applications in medicine, such as drug delivery, lithotripsy, hyperthermia and diagnostic imaging [ 5 - 10 ]. Ultrasonic medical applications lead to improving treatment efficiency and decreasing side effects. In addition, the ultrasonic devices help to achieve more information in a non-invasive diagnosis process.
There are two main approaches for ultrasonic lithotripsy; Shock Wave Lithotripsy (SWL) machines typically use single-cycle pulses at a low frequency and high peak pressures. In contrast, Burst Wave Lithotripsy (BWL) applies short bursts of focused, sinusoidal ultrasonic pulses [ 11 ]. Shock wave lithotripsy remains a popular treatment and is the only non-invasive lithotripsy in clinical application since 1980 [ 12 , 13 ]. The SWL, as an effective biomedical treatment, is the most common treatment for kidney stones and also a useful method for uncomplicated, upper urinary tract calculi. The shock wave’s effects and cavitation phenomenon are two important mechanisms on kidney stone fragmentation [ 14 ]. K. G. Wang [ 15 ] explained the effects of a bubble in the shock path and fracture behaviors of the stone using a multiphase fluid-solid coupled model. Songlin Zhu et al. [ 16 ] simulated the stone fragmentation in the renal pelvis and investigated the role of stress waves and cavitation in stone comminution for shock-wave lithotripsy (SWL).
Several researches have been done on lithotripsy of kidney stones. Chaussy et al. [ 17 ] investigated the fracture of calculi using focused shock waves instead of surgery and reduced the need for surgery with the aid of shock wave lithotripsy. Xi and zhong [ 18 ] investigated the transient stress fields, produced in stones with different geometry and size, and the general pattern of wave propagation using a dynamic photoelastic imaging technique. Ying Zhang et al. [ 19 ] studied the effect of stone size on the comminution process and efficiency in SWL. Dahake and Gracewski [ 20 , 21 ] developed a linear elastic model of stress waves within the stone and studied the spherical and oval models of stone. Cleveland and Sapozhnikov [ 22 ] applied a linear elasticity model based on the finite-difference method to study the fragmentation of kidney stone, and determined the location of maximum stress. Weinberg and Eritz [ 23 ] numerically studied the fracture of kidney stones. They compared two shock waves with the same energy transport and found that the wave with larger tensile amplitude has more fracture influence.
Most numerical studies on lithotripsy have been done by applying the elasticity method or energy method and neglected the dissipation phenomenon. In this research, the effects of ultrasonic waves on the stone and surrounding tissues are investigated using the non-linear Westervelt equation. First, the comparison between the application of linear and nonlinear waves is provided. Then, the von Mises stress and deformation in stone, the effect of stone material and thermal analysis of SWL on the stone and surrounding tissues are studied. Thermal analysis asserts information about side effects and damage on surrounding tissue.
Material and Methods
This computational study is done by COMSOL Multiphysics. It is based on the ﬁnite element method for discretization and solution of governing. The nature of the problem requires the combination and coupling of acoustic wave propagation, elasticity of structure and biological heat transfer equations. One-way coupling is employed to couple the acoustics and structural physics.
The schematic of the proposed model, kidney stone and surrounding tissues is shown in Figure 1. Surrounding tissues, including skin, fat layer, kidney tissue and water, are considered in the simulation. Space around the stone is filled with urine [ 10 ]. Urine is composed of 95% water [ 24 ]. Therefore, water is used for the surroundings of the stone. Tissues properties are provided in Table 1.
|Tissue||Speed of Sound c (m/s)||Density ρ(kg/m3)||Nonlinearity Parameter B/A||Sound Absorption Coefficient α (Np/m/MHz)||Specific Heat Capacity C (J/Kg.c°)||Thermal Conductivity k (W/m.c°)|
|COM: Calcium Oxalate Monohydrate|
In all simulations, the material of stone is Calcium Oxalate Monohydrate (COM). In addition, Uric acid is used to specify the effect of stone material. The Stones properties are given in Table 2.
|Stone||Density (kg/m3)||Speed of Sound (m/s)||Young's Modu-lus (GPa)||Poisson's Ratio||Tensile Strength (MPa)|
|Calcium Oxalate Monohydrate (COM)||2038||4535||24.51||0.333||1.1|
Investigation of SWL lithotripsy and its thermal side effects on surrounding tissues are the Multiphysics problem and need to couple acoustics, elasticity and bio-heat transfer equations. The governing equations are as follows.
The linear elastic wave equation (LEWE) model describes wave propagation within three dimensional, linear, isotropic elastic media. The former is formulated as shown below:
ρ ∂ 2 v x ∂ t 2 = ∂ τ xx ∂ x + ∂ τ xy ∂ y + ∂ τ xz ∂ z + f x (1)
ρ ∂ 2 v y ∂ t 2 = ∂ τ xy ∂ x + ∂ τ yy ∂ y + ∂ τ yz ∂ z + f y (2)
ρ ∂ 2 v z ∂ t 2 = ∂ τ zy ∂ x + ∂ τ zy ∂ y + ∂ τ zz ∂ z + f z (3)
where ρ is the medium density, vx, vy, vz are the displacement components along the x, y and z directions, respectively, ∂τxx, ∂τyy, ∂τzz and ∂τxy, ∂τxz, ∂τyz are the normal and shear stress components, respectively, and fx, fy, fz are the body-force components. The stress-strain relations are:
τ xx = ( λ + 2 μ ) ∂ v x ∂ x + λ ( ∂ v y ∂ y + ∂ v z ∂ z ) (4)
τ yy = ( λ + 2 μ ) ∂ v y ∂ y + λ ( ∂ v x ∂ x + ∂ v z ∂ z ) (5)
τ zz = ( λ + 2 μ ) ∂ v z ∂ z + λ ( ∂ v x ∂ x + ∂ v y ∂ y ) (6)
τ xy = μ ( ∂ v x ∂ y + ∂ v y ∂ x ) (7)
τ xz = μ ( ∂ v x ∂ z + ∂ v z ∂ x ) (8)
τ yz = μ ( ∂ v y ∂ z + ∂ v z ∂ y ) (9)
Where λ and μ are the Lam´e coefficients.
The Westervelt-Lighthill Equation (WLE) is a nonlinear full-wave equation mostly used to model acoustic field propagation in nonlinear thermoviscous fluids. The Westervelt PDE is as follows:
∇ 2 p - 1 c 2 ∂ 2 p ∂ t 2 + δ c 4 ∂ 3 p ∂ t 3 = - β ρ 0 c 4 ∂ 2 P 2 ∂ t 2 (10)
Where p is the acoustic pressure, c is the speed of sound and ρ0 is the ambient density. In addition, δ and β are the sound diffusivity of medium and nonlinearity coefficient, respectively, which are defined in Eq. 11
β = 1 + B 2 A , δ = 2 α c 3 ω 3 (11).
The linear elasticity equation is applied.
ρ ∂ 2 u ∂ t 3 = ∇s + F V (12)
Where u is displacement and FV is external force, solved and added by the acoustic equation.
The Pennes Bioheat Equation is used to calculate the thermal changes and temperature distribution in tissues.
ρ 0 C t ∂T ∂t = ∇k ∇T + C b w b ( T - T b ) + Q abs (13)
Q abs = 2 α abs I = α P 2 ρc (14)
After solving the acoustic equation in tissue, heat generation of the ultrasonic waves (acoustic intensity) can be calculated.
In the Eq. 13, k is the thermal conductivity of tissue; wb is the flow rate of blood through the tissue; Cb is the specific heat capacity of blood. Qabs is acoustic intensity, which is equivalent to the absorbed ultrasonic power in the tissue. Ct and ρ are specific heat capacity and density of the tissue, respectively. The perfusion time in the blood is in order of several thousand seconds and can be neglected for short period.
In Eq. 14, α is the absorption coefficient; P and ρ are acoustic pressure and tissue density, respectively.
The pressure inlet is applied for the transducer boundary. Upper, right and left walls are considered as a perfectly matched layer.
Background (initial) pressure within the entire domain is set to zero and the initial temperature is set to 37 °C.
In acoustic simulations, considering the order of discretion, the number of meshes usually should be chosen at least 4 meshes per wavelength. For the correct resolution in the focal region and the places where the pressure gradient is higher such as kidney and stones, the maximum mesh size is set to λ/6 with discrete order of 4 and in other regions λ/4 with discrete order of 2. The number of meshes depends on the maximum wavelength of the shock wave.
Mesh independency is studied for Church impulse profile by three different sizes of mesh. Trian-gular meshes are used, which the numbers of elements are 69958, 76887 and 83354, respectively. The pressure variety versus time at the front of the stone for different meshes is shown in Figure 2.
As shown in Figure 2, the results for 2nd and 3rd mesh have a good match. Therefore, the 2nd mesh is applied for simulation.
Wienberg and Ortiz [ 23 ] simulated SWL on a domain with constant properties and measured the pressure through it. Their results are used for validation. Their proposed model is depicted in Figure 3a.
To calibrate the simulation, the shock wave is emitted in a homogeneous elastic media with 1 MPa of Young’s modulus and 1540 m/s sound velocity. The shock values are set to τ1=1.1 μs, τ2=1.96 μs, P=100 MPa and the wave diagram is measured at 10 and 40 mm.
Wienberg and Ortiz used the energy method and applied Eq. 15 to consider the dissipation parameter but we applied the Westervelt acoustic method with β=4.7, α = N P m * MHz .
P damp = α ρ c L l e ϑ . (15)
Where α is the damping coefficient; ρ is density, CL is the speed of the sound; le is the longitudinal character of the element and ϑ . is volume traction rate.
The comparison of the current study and the Wienberg’s results are presented in Figure 3b.
As depicted in Figure 3b, the overall trend and advent time of waves are similar. The difference of dissipation consideration leads to a difference in maximum pressure of the 40 mm data series between the current study and the Wienberg’s results. At high pressures and frequencies, it is expected that wave misshape its initial form due to its nonlinearity, and distorts to a saw blade shape. This point is clearly specified by comparisons of 40 mm and 10 mm series in our simulation. In contrast, it is not indicated in the Wienberg’s results. Therefore, the results of the current study are more accurate.
Comparison between Westervelt and elastic (linear) methods
The maximum pressure in the linear method, which the sound absorption is not considered, is about 91 MPa, and this parameter in the Westervelt method is 28.6 MPa. The wave arrival time to the stone for both simulations is 28.5 microseconds. In addition, the transverse wave amplitudes in the linear solution and the Westervelt solution are about 10 mm and 20 mm, respectively.
The acoustic pressure when the wave reaches and contacts with stone for both methods are shown in Figure 4.
Pressure profiles of linear and Westervelt method along z-axis are depicted in Figure 5.
As shown in Figure 5, the highest pressure at the moment of the wave reaching the stone in linear solution is at Z= 42.2 mm and in a Westervelt solution is at Z= 45 mm; this difference is due to the consideration of nonlinearity parameter in the Westervelt method. As mentioned, because of the consideration of the dissipation phenomenon, the highest pressure in the Westervelt method is minor.
The pressure is maximum at the front of the stone, where the velocity decreases due to properties changing. High density and high speed of sound cause to accumulate the pressure at the front of the stone.
Distribution of von Mises stress and stone deformation for both methods are depicted in Figure 6.
Figure 6 shows the stress and deformation for the highest pressure in the stone. The maximum pressure occurs at t=35 µs, which is the same in both linear and Westervelt methods.
Distributions of acoustic intensity for both methods are depicted in Figure 7.
The acoustic wave leads to different levels of injury in tissue. The destructive thresholds of energy flux (acoustic intensity) on the tissues are given in Table 3.
|Energy Flux (mJ/mm2)||Biologic Effects|
|0.3||The breakdown of endothelial cell structure usually occurs at the inlet of pulse pressure due to positive pressure|
|0.22||Changes in subcellular structures especially mitochondria|
|0.1||Formation of tensile structure in the endothelium|
|0.007||Damage to Skin|
According to Figure 7 and Table 3, the shock wave damages surrounding tissues in both linear and nonlinear simulation, but acoustic intensity significantly is greater in the linear method than the nonlinear method. The damaged area, obtained in our results, is in the focal region and before it, which has good matching with experimental and medical consequences [ 23 ].
Effect of shock profiles
To compare different impulse profiles of shock waves and study the effect of the rise and inactivation times on the produced acoustic pressure, pressure amplitude is set to 35 MPa. Different impulse profiles are shown in Figure 8.
As depicted in Figure 8, Church impulse and Cleveland & Sapozhnikov impulse have similar rise times. The profile of Bailey impulse approximately is similar to Colonius impulse.
Acoustics waves propagations through the stone and tissues are shown in Figure 9 for Bailey, Church and Energy equivalent Church impulses.
According to Figures 8 and 9, the acoustic pressure increases by decreasing the rise time. Moreover, by reduction of the rise time, the focal zone becomes larger and as a result, the fracture occurs in more parts of the stone. However, enlargement of focal zone damages to the tissues and causes hematoma in the kidney.
In this section, Uric acid is utilized as well as Calcium Oxalate Monohydrate. These two stones are the most common kidney stones. Von Mises stress in both stones is depicted in Figure 10.
As shown in Figure 10, the maximum pressure in Uric acid stone is smaller than Calcium Oxalate Monohydrate but the maximum pressure region is more extensive; thus, it can be expected with the same acoustics pressure, the fracture of Uric Acid occurs in a larger area.
The thermal simulation of a shock wave is done for 64 microseconds. The time of shock wave radiation is very short, about 10-20 microseconds and the frequency of repetition is 1-2 times per second. The maximum produced power per shock is more than 1 MW. The average energy per pulse is 10-150 mJ in the focal region. The temperature distribution is shown in Figure 11.
As depicted in Figure 11, thermal analysis shows that the temperature about a few thousandths degrees has risen in stone and no effect on the tissues.
In the current study, the non-linear Westervelt equation is used to simulate the ultrasonic shockwave for lithotripsy purposes on Calcium Oxalate Monohydrate stone and surrounding tissues. In contrast, most numerical studies in lithotripsy have been carried out using the elasticity or energy methods and neglected the dissipation phenomenon. First, the comparison between linear and nonlinear waves is provided. As shown, the highest pressure decreases by the involvement of the dissipation phenomenon. In fact, the highest pressure in the Westervelt method is less than the elastic method. By comparing acoustic intensity with destructive thresholds of energy flux on the tissues for both methods, it is indicated that shock wave injures kidney tissue. The damage and its location have a good match with experimental and medical consequences [ 23 ]. Then, the different impulse profiles of shock waves are compared to specify the effect of the rise and inactivation times on the produced pressure. As depicted, the produced pressure increases by decreasing the rise time. Next, Uric acid is applied as well as Calcium Oxalate Monohydrate to explain the effect of stones’ material on the stress and its distribution. There is a significant difference, which asserts the stone’s material has a key role in fragmentation. Finally, the temperature distribution is obtained by the bioheat equation. As illustrated, the rise in temperature is minimal in shock wave lithotripsy and the resulting damage is very negligible compared to the HIFU. The difference of emitting time leads to this variety of damage; thus, that emitting time is about a few microseconds in fragmentation of kidney stone but it is a few seconds in the HIFU.
In this study, the ultrasonic shock waves Lithotripsy and its effects on surrounding tissues are investigated by coupling and solving the acoustic wave propagation, elasticity of structure and biological heat transfer equations.
To conclude, the following results are obtained:
- • The effect of linear and nonlinear wave consideration is discussed and found that disaffiliation of the wave nonlinearity causes a high incompatibility with reality.
- • The shock wave damages surrounding tissues in both linear and nonlinear simulation, but acoustic intensity significantly is greater in the linear method than the nonlinear method. The damaged area is the focal region and before it.
- • By comparing different impulse profiles of shock waves, it is stated the produced acoustic pressure increases by decreasing the rise time. Moreover, by reduction of the rise time, the focal zone becomes larger and as a result, the fracture occurs in more parts of the stone.
- • By changing the stone material, it is declared with the same acoustic pressure, the fracture of Uric Acid occurs in a larger area than Calcium Oxalate Monohydrate.
- • Thermal analysis of shock wave lithotripsy shows that the temperature variation is minimal and the resulting damage of incremental temperature is very negligible compared to the HIFU.
Azadeh Shahidian conceived the idea. All authors contributed to the study conception and design. The method implementation was carried out by Mahdi Moghim Nejad and Azadeh Shahidian. Material preparation, data collection and analysis were performed by Mahdi Moghim Nejad and Azadeh Shahidian. Also, the first draft of the manuscript was written by them and all authors commented on previous versions of the manuscript. All authors read and approved the final version of the manuscript.
Conflict of Interest
- Mohammadalibeigi F, Shirani M, Seyed-Salehi H, Afzali L. Biochemical urinalysis of healthy kidney and stone-generating kidney in unilateral urolithiasis. Journal of Renal Injury Prevention. 2019; 8(2):151-6. DOI
- Zeng J, Wang S, Zhong L, Huang Z, Zeng Y, Zheng D, Zou W, Lai H. A Retrospective Study of Kidney Stone Recurrence in Adults. Journal of Clinical Medicine Research. 2019; 11(3):208-12. Publisher Full Text | DOI | PubMed
- Eliahou R, Hidas G, Duvdevani M, Sosna J. Determination of renal stone composition with dual-energy computed tomography: an emerging application. Seminars in Ultrasound, CT and MRI. 2010; 31(4):315-20. DOI
- Warty Y, Haryanto F, Fitri LA, Haekal M, Herman H. A Spatial Distribution Analysis on the Deposition Mechanism Complexity of the Organic Material of Kidney Stone. J Biomed Phys Eng. 2020; 10(3):273-82. DOI
- Ueberle F. Application of shock waves and pressure pulses in medicine. Springer Handbook of Medical Technology: Berlin, Heidelberg; 2011.
- Beik J, Mehdizadeh AR, Shakeri-Zadeh A. Ultrasound in cancer treatment through nanotechnology. J Biomed Phys Eng. 2016; 6(3):123-6. Publisher Full Text | PubMed
- Ghassemi M, Shahidian A. Nano and bio heat transfer and fluid flow. 2017.
- Aayani R, Shahidian A, Ghassemi M. Numerical Investigation of Non-Newtonian Blood Effect on Acoustic Streaming. Journal of Applied Fluid Mechanics. 2016; 9(1):173-6.
- López-Haro SA, Gutiérrez MI, Vera A, Leija L. Modeling the thermo-acoustic effects of thermal-dependent speed of sound and acoustic absorption of biological tissues during focused ultrasound hyperthermia. Journal of Medical Ultrasonics. 2015; 42(4):489-98. DOI
- Bailey MR, Khokhlova VA, Sapozhnikov OA, Kargl SG, Crum LA. Physical mechanisms of the therapeutic effect of ultrasound (a review). Acoustical Physics. 2003; 49(4):369-88. DOI
- Dai JC, Bailey MR, Sorensen MD, Harper JD. Innovations in Ultrasound Technology in the Management of Kidney Stones. Urologic Clinics. 2019; 46(2):273-85. DOI
- Wess OJ, Mayer J. Fragmentation of brittle material by shock wave lithotripsy. Momentum transfer and inertia: a novel view on fragmentation mechanisms. Urolithiasis. 2018;137-49. DOI
- Lawler AC, Ghiraldi EM, Tong C, Friedlander JI. Extracorporeal shock wave therapy: current perspectives and future directions. Current Urology reports. 2017; 18(4):25. DOI
- Ghorbani M, Oral O, Ekici S, Gozuacik D, Koşar A. Review on lithotripsy and cavitation in urinary stone therapy. IEEE Reviews in Biomedical Engineering. 2016; 9:264-83. DOI
- Wang KG. Multiphase fluid-solid coupled analysis of shock-bubble-stone interaction in shockwave lithotripsy. International Journal for Numerical Methods in Biomedical Engineering. 2017; 33(10):e2855. DOI
- Zhu S, Cocks FH, Preminger GM, Zhong P. The role of stress waves and cavitation in stone comminution in shock wave lithotripsy. Ultrasound in Medicine & Biology. 2002; 28(5):661-71. DOI
- Chaussy CH, Brendel W, Schmiedt E. Extracorporeally induced destruction of kidney stones by shock waves. The Lancet. 1980; 316(8207):1265-8. DOI
- Xi X, Zhong P. Dynamic photoelastic study of the transient stress field in solids during shock wave lithotripsy. J Acoust Soc Am. 2001; 109(3):1226-39. DOI
- Zhang Y, Nault I, Mitran S, Iversen ES, Zhong P. Effects of stone size on the comminution process and efficiency in shock wave lithotripsy. Ultrasound in Medicine & Biology. 2016; 42(11):2662-75. DOI
- Dahake G, Gracewski SM. Finite difference predictions of P-SV wave propagation inside submerged solids. I. Liquid–solid interface conditions. J Acoust Soc Am. 1997; 102(4):2125-37. DOI
- Dahake G, Gracewski SM. Finite difference predictions of P-SV wave propagation inside submerged solids. II. Effect of geometry. J Acoust Soc Am. 1997; 102(4):2138-45. DOI
- Cleveland RO, Sapozhnikov OA. Modeling elastic wave propagation in kidney stones with application to shock wave lithotripsy. J Acoust Soc Am. 2005; 118(4):2667-76. DOI
- Weinberg K, Ortiz M. Kidney damage in extracorporeal shock wave lithotripsy: a numerical approach for different shock profiles. Biomechanics and Modeling in Mechanobiology. 2009; 8(4):285. DOI
- Handbook of LC-MS bioanalysis: best practices, experimental protocols, and regulations. 2013.
- Kyriakou A. Multi-physics computational modeling of focused ultrasound therapies. 2015.
- Suomi V, Jaros J, Treeby B, Cleveland RO. Full Modeling of High-Intensity Focused Ultrasound and Thermal Heating in the Kidney Using Realistic Patient Models. IEEE Trans Biomed Eng. 2018; 65(11):2660-70. DOI
- Haddadi S, Ahmadian MT. Numerical and Experimental Evaluation of High-Intensity Focused Ultrasound–Induced Lesions in Liver Tissue Ex Vivo. Journal of Ultrasound in Medicine. 2018; 37(6):1481-91. DOI
- Nyame YA, De S, Sarkissian C, Brown R, Kartha G, Babbar P, Monga M. Kidney stone models for in vitro lithotripsy research: a comprehensive review. Journal of Endourology. 2015; 29(10):1106-9. DOI
- Steinbach P, Wörle K, Seidl M, Seitz R, Hofstädter F. Effects of high-energy ultrasonic shock waves on tumor cells in vitro and human endothelial cells in situ shock wave lithotripsy, aspects and Forecasts. University Library: Tübingen; 1993.
- Miller DL, Thomas RM. Thresholds for hemorrhages in mouse skin and intestine induced by lithotripter shock waves. Ultrasound in Medicine & Biology. 1995; 21(2):249-57. DOI