k-Wave
A MATLAB toolbox for the time-domain
simulation of acoustic wave fields
- Getting Started
- Examples
- Initial Value Problems
- Example: Homogenous Propagation Medium
- Example: Using A Binary Sensor Mask
- Example: Defining A Sensor Mask By Opposing Corners
- Example: Loading External Image Maps
- Example: Heterogeneous Propagation Medium
- Example: Saving Movie Files
- Example: Recording The Particle Velocity
- Example: Defining A Gaussian Sensor Frequency Response
- Example: Comparison Of Modelling Functions
- Example: Setting An Initial Pressure Gradient
- Example: Simulations In One Dimension
- Example: Simulations In Three Dimensions
- Example: Photoacoustic Waveforms in 1D, 2D and 3D
- Time Varying Source Problems
- Example: Monopole Point Source In A Homogeneous Propagation Medium
- Example: Dipole Point Source In A Homogeneous Propagation Medium
- Example: Simulating Transducer Field Patterns
- Example: Steering A Linear Array
- Example: Snell's Law And Critical Angle Reflection
- Example: The Doppler Effect
- Example: Diffraction Through A Slit
- Example: Simulations In Three-Dimensions
- Sensor Directivity
- Example: Focussed Detector in 2D
- Example: Focussed Detector in 3D
- Example: Modelling Sensor Directivity in 2D
- Example: Modelling Sensor Directivity in 3D
- Example: Sensor Element Directivity in 2D
- Example: Focussed 2D Array with Directional Elements
- Photoacoustic Image Reconstruction
- Example: 2D FFT Reconstruction For A Line Sensor
- Example: 3D FFT Reconstruction For A Planar Sensor
- Example: 2D Time Reversal For A Line Sensor
- Example: 2D Time Reversal For A Circular Sensor
- Example: 3D Time Reversal For A Planar Sensor
- Example: 3D Time Reversal For A Spherical Sensor
- Example: Image Reconstruction With Directional Sensors
- Example: Image Reconstruction With Bandlimited Sensors
- Example: Iterative Image Improvement Using Time Reversal
- Example: Attenuation Compensation Using Time Reversal
- Example: Attenuation Compensation Using Time Variant Filtering
- Example: Automatic Sound Speed Selection
- Diagnostic Ultrasound Simulation
- Example: Defining An Ultrasound Transducer
- Example: Simulating Ultrasound Beam Patterns
- Example: Using An Ultrasound Transducer As A Sensor
- Example: Simulating B-mode Ultrasound Images
- Example: Simulating B-mode Images Using A Phased Array
- Numerical Analysis
- Example: Controlling The Absorbing Boundary Layer
- Example: Source Smoothing
- Example: Filtering A Delta Function Input Signal
- Example: Modelling Power Law Absorption
- Example: Modelling Nonlinear Wave Propagation
- Example: Optimising k-Wave Performance
- Using The C++ Code
- Elastic Wave Propagation
- Example: Explosive Source In A Layered Medium
- Example: Plane Wave Absorption
- Example: Shear Waves And Critical Angle Reflection
- Example: Simulations In Three Dimensions
- Functions - By Category
- Functions - Alphabetical List
- Release Notes
- License
k-Wave Toolbox |
Controlling The Absorbing Boundary Layer Example
On this page… |
---|
Overview
The first-order simulation codes included within the k-Wave Toolbox (kspaceFirstOrder1D
, kspaceFirstOrder2D
, and kspaceFirstOrder3D
) use a special type of anisotropic absorbing boundary layer known as a perfectly matched layer (PML) to absorb acoustic waves when they reach the edge of the computational domain. By default, this layer occupies a strip of 20 grid points (10 in 3D) around the edge of the domain. Without this boundary layer, the computation of the spatial derivates via the FFT causes waves leaving one side of the domain to reappear at the opposite side. The use of the PML thus facilitates infinite domain simulations without the need to increase in the size of the computational grid. This example demonstrates how to control the parameters of the PML within k-Wave via optional input parameters.
Controlling the properties of the absorbing boundary layer
The PML has four properties, each of which can be controlled via optional input parameters:
- Size: The size of the layer on each edge of the domain is set by
'PMLSize'
in units of grid points. The default setting is 20 grid points in 1D and 2D and 10 grid points in 3D. If the size is specified as a single number, this is used for all Cartesian directions. Alternatively, the size for each direction can be set individually by setting'PMLSize'
to[PML_X_SIZE, PML_Y_SIZE]
in 2D, and[PML_X_SIZE, PML_Y_SIZE, PML_Z_SIZE]
in 3D. - Absorption: The absorption within the layer is set by
'PMLAlpha'
in units of Nepers per grid point. The default setting is 2 for all dimensions. Again, if the absorption is specified as a single number, this is used for all Cartesian directions. Alternatively, the absorption for each direction can be set individually by setting'PMLAlpha'
to[PML_X_ALPHA, PML_Y_ALPHA]
in 2D, and[PML_X_ALPHA, PML_Y_ALPHA, PML_Z_ALPHA]
in 3D. - Location: The location of the layer can be set to be inside or outside the computational grid created by the user by changing the value of the Boolean flag
'PMLInside'
(the default setting istrue
). If'PMLInside'
is set tofalse
, the input grids are enlarged on each edge by the size given by'PMLSize'
. - Visibility: The visibility of the absorbing layer within the animations displayed during the simulation is controlled by the Boolean input parameter
'PlotPML'
(the default setting istrue
).
For accurate simulations, care must be taken that the initial pressure distribution and the sensor mask don't lie within the PML. This can be avoided by setting 'PMLInside'
to false
. However, the computational time will still be dependent on the total size of the grid including the PML (i.e., computations will be fastest for grids where the total number of grid points in each dimension is given by a power of two).
Switching off the PML
The effect of the PML can be illustrated most clearly by switching it completely off by setting the optional input 'PMLAlpha'
to 0
(set example_number = 1
within the example m-file). This sets the absorption within the layer to zero, thus the waves leaving one side of the domain reappear at the opposite side. A visualisation of the recorded time-series is given below.
A similar effect is seen if the value for 'PMLAlpha'
is set too high (set example_number = 2
within the example m-file). This causes the waves to be reflected from the layer.
A partially effective PML
The effectiveness of the PML depends on its size and absorption, but also on the time step used in the simulation (the more time steps the wave spends in the layer, the more it will be absorbed). A partially effective PML can be illustrated by reducing the size of the PML by setting the optional input 'PMLSize'
to 2
(set example_number = 3
within the example m-file). The layer absorbs some (but not all) of the waves approaching the boundaries, and some of the wave wrapping seen in the previous example is still visible.
Setting the PML to be outside the computational grid
By default, the PML will be inside the grid defined by the user, so care must be taken that the source and sensor are not defined within this region. This can be avoided by setting the optional input 'PMLInside'
to false
(set example_number = 4
within the example m-file). In this case the computational grid and input variables are automatically expanded by the size of the PML. This setting gives several additional command line status readouts, including feedback on expanding the input grids and the total size of the computational grid (in this case 168 by 168 grid points).
Running k-Wave simulation... start time: 02-Nov-2011 20:16:41 reference sound speed: 1500m/s dt: 20ns, t_end: 12.06us, time steps: 604 input grid size: 128 by 128 grid points (12.8 by 12.8mm) maximum supported frequency: 7.5MHz smoothing p0 distribution... expanding computational grid... computational grid size: 168 by 168 grid points calculating Delaunay triangulation... precomputation completed in 1.1535s starting time loop... estimated simulation time 5.9372s... simulation completed in 9.6409s total computation time 10.8585s
As the inner part of the computation grid is identical to that in the Homogeneous Propagation Medium Example, the recorded time-series (shown below) is also identical. In this case, the default size and absorption for the PML are adequate to prevent any visible wave wrapping. For any particular example, if the performance of the PML is not sufficient, increasing the size of the PML, optimising the value of PMLAlpha
, or reducing the time-step used in t_array
will give improvements.
Numerical Analysis | Source Smoothing |
© 2009-2014 Bradley Treeby and Ben Cox.