Full-Wave Modeling of Transcranial Ultrasound using Volume-Surface Integral Equations and CT-Derived Heterogeneous Skull Data

TL;DR

This study develops a velocity-stress–based full-wave solver using volume-surface integral equations (VSIE) to simulate transcranial ultrasound through heterogeneous skull bone directly from CT voxel data. By comparing VSIE with high-accuracy boundary element method (BEM) results and performing mesh-convergence studies, the authors demonstrate accurate focal predictions even with voxel-sized meshes, and they quantify how skull heterogeneity distorts and attenuates the focus. Key findings show substantial focal aberration and elevated intra-bone pressures when skull heterogeneity is included, as well as good agreement between VSIE and BEM in homogeneous scenarios, validating the approach. The methodology offers a data-preserving pipeline from CT to acoustic parameters and supports six-voxel-per-wavelength accuracy, with clear implications for safer and more effective transcranial ultrasound planning.

Abstract

Transcranial ultrasound therapy uses focused acoustic energy to induce therapeutic bioeffects in the brain. Ultrasound must be transmitted through the skull, which is highly attenuating and heterogeneous, causing beam distortion, reducing focal pressure, and shifting the target location. Computational models are frequently used to predict beam aberration, assess cranial heating, and correct the phase of ultrasound transducers. These models often rely on computed tomography (CT) images to build patient-specific geometries and estimate skull acoustic properties. However, the coarse voxel resolution of CT limits accuracy for differential equation solvers at ultrasound frequencies. This paper presents an efficient numerical method based on volume-surface integral equations to model full-wave acoustic propagation through heterogeneous skull bone. We show that our approach effectively simulates transcranial ultrasound, even when using the original CT voxels as the computational mesh, where the 0.5 mm voxel length is relatively coarse compared to the shortest wavelength of 3 mm. The method is validated against a high-resolution boundary element model using an averaged skull representation. Simulations using a CT-based skull model and a bowl transducer reveal significant beam distortion of 7.8 mm attributed to the skull's heterogeneous acoustical properties.

Figures (9)

• Figure 1: A sketch of the simulation domain, which consists of a transducer source $s$, a bone region $\Omega_\mathrm{bone}$ with normal $\hat{\mathbf{n}}$ at its surface $\Gamma$, and an unbounded exterior tissue domain $\Omega_\mathrm{ext}$.
• Figure 2: A visual interpretation of the algorithmic pipeline to calculate density and sound speed values inside skull bone from CT images.
• Figure 3: The histograms and boxplots show the distribution of the voxel values concerning the density and speed of sound. The computational grid of the skull slab consists of 199,693 voxels.
• Figure 4: The geometry of the ultrasound transducer and the skull slab.
• Figure 5: The magnitude of the acoustic field for the four cases: the incident field emitted by the bowl transducer (top row), the total field simulated with the VSIE and heterogeneous bone (second row), the VSIE with homogeneous bone (third row), and the BEM simulation with homogeneous bone (bottom row). The pressure field is calculated on different slices on the $x$-axis. The center of the bowl transducer is located at the global origin and emits an acoustic field of 500 kHz in the negative $y$-direction.