Example 1ΒΆ

The example shows how to prepare k-Wave input file describing simulation in linear homogeneous medium. The pressure source is specified in the form of Initial pressure distribution (\(p_0\)).

The simulation has spatial resolution of \(10^{-3}\) m and the simulation domain of \(64 \times 128 \times 256\) grid points. While temporal resolution is computed so that CFL condition number is 0.3.

Finally, the sensor mask is set up to allow collection of pressure evolution in time at a plane perpendicular to Z axis placed in the middle of the domain (in Z-direction).

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
import sys
import numpy as np

sys.path.append('..')

from kwave_input_file import KWaveInputFile
from kwave_data_filters import SpectralDataFilter

# Define input file name, domain resolution and simulation grid size
filename = 'example_001_simple_p0_py.h5'
grid_res = (1.0e-4, 1.0e-4, 1.0e-4)
grid_size = (64, 128, 256)

# Compute simulation time-step length, choose simulation length
# and compute number of time-steps required
dt = 0.3 * grid_res[0] / 1500.0
end_time = ((grid_size[0] / 2) * grid_res[0]) / 1500.0
steps = int(end_time / dt)

# Open k-Wave input file object
with KWaveInputFile(filename, grid_size, steps, grid_res, dt) as file:
    # Write homogeneous medium sound-speed of 1500 m/s into the file
    file.write_medium_sound_speed(1500.0)
    # Write homogeneous medium density of 1000 kg/m^3 into the file
    file.write_medium_density(1000.0)

    # Prepare initial pressure distribution in the domain
    # (16^3 gp cube of 1kPa in the middle of the domain)
    p0 = np.zeros(grid_size)
    p0[int(grid_size[0] / 2 - 16): int(grid_size[0] / 2 + 16),
       int(grid_size[1] / 2 - 16): int(grid_size[1] / 2 + 16),
       int(grid_size[2] / 2 - 16): int(grid_size[2] / 2 + 16)] = 1.0e3
    # Smooth the initial pressure using Blackman window in spectral filter
    p0 = SpectralDataFilter.smooth(p0)
    # Write the initial pressure into the file
    file.write_source_input_p0(p0)

    # Setup pressure sampling mask as X-Y aligned plane in the middle of
    # the domain
    sensor_mask = np.zeros(grid_size)
    sensor_mask[:, :, int(grid_size[2] / 2)] = 1
    # Transform the mask into indexes in the ordering of the output data
    file.write_sensor_mask_index(file.domain_mask_to_index(sensor_mask))