k-Wave User Forum » Topic: linear array UT simulation in heterogeneous medium

k-Wave User Forum » Topic: linear array UT simulation in heterogeneous medium http://www.k-wave.org/forum/topic/linear-array-ut-simulation-in-heterogeneous-medium Support for the k-Wave MATLAB toolbox en-US Wed, 13 May 2026 00:26:13 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "linear array UT simulation in heterogeneous medium" http://www.k-wave.org/forum/topic/linear-array-ut-simulation-in-heterogeneous-medium#post-4136 Sun, 24 Nov 2013 15:20:51 +0000 Bradley Treeby 4136@http://www.k-wave.org/forum/ Hi Eun, There seems to be a couple of issues. First, your perfectly matched layer (PML) is too small to be effective. This should be at least 10 grid points, otherwise you will notice the waves wrapping from one side of the domain to the other. There's some more information in Section 2.6 of the <a href="http://www.k-wave.org/documentation.php">k-Wave Manual</a> that might be helpful (see also Fig. 2.3), as well as in the Numerical Analysis -> Controlling The Absorbing Boundary Layer Example in the MATLAB documentation. Second, it appears that your transducer is within the PML. You should make sure that your source is outside the PML, or set <code>'PMLInside'</code> to <code>false</code>. Remember, your simulations will be fastest when the total grid dimensions (including the PML) have small prime factors. You can check this using <code>checkFactors</code>. Brad. sheun on "linear array UT simulation in heterogeneous medium" http://www.k-wave.org/forum/topic/linear-array-ut-simulation-in-heterogeneous-medium#post-4098 Fri, 22 Nov 2013 02:37:49 +0000 sheun 4098@http://www.k-wave.org/forum/ first i would like to thank you for the tool, it helped a lot in understanding acoustics in general. I am trying to model linear array sensor in heterogenous medium and see how the reflected signal from the medium boundary looks like. Two mediums with different "sound speed", "density", "att. coefficient" are located top and bottom, and the sensor are located at the top surface of the upper medium. When I run the code, it looks like the waves are propagating in both direction. I was expecting for the waves to propagate in only one direction (propagating towards the bottom). This doesn't happen when I switch the medium parameters between the top and bottom. Below is the code i used. I really appreciate if you could give me some direction as to how I should modify the codes. clear all; % simulation settings DATA_CAST = 'single'; % ========================================================================= % DEFINE THE K-WAVE GRID % ========================================================================= % set the size of the perfectly matched layer (PML) PML_X_SIZE = 2; % [grid points] PML_Y_SIZE = 1; % [grid points] PML_Z_SIZE = 1; % [grid points] % set total number of grid points not including the PML Nx = 196 - 2*PML_X_SIZE; % [grid points] Ny = 126 - 2*PML_Y_SIZE; % [grid points] Nz = 62 - 2*PML_Z_SIZE; % [grid points] % set desired grid size in the x-direction not including the PML x = 0.75/1000; % [m] % calculate the spacing between the grid points dx = x; % [m] dy = dx; % [m] dz = dx; % [m] % create the k-space grid kgrid = makeGrid(Nx, dx, Ny, dy, Nz, dz); % ========================================================================= % DEFINE THE MEDIUM PARAMETERS % ========================================================================= % define the properties of the propagation medium medium.sound_speed = 2000*ones(Nx,Ny,Nz) ; % [m/s] 2nd medium.sound_speed(1:Nx/2,:,:) = 1000 ; % [m/s] 1st 101 medium.density = 2000*ones(Nx,Ny,Nz); % [kg/m^3] 2nd medium.density(1:round(Nx/2),:,:) = 1134 ; % [kg/m^3] 1st medium.alpha_coeff = 0.001*ones(Nx,Ny,Nz); % [dB/(MHz^y cm)] 2nd medium.alpha_coeff(1:round(Nx/2),:,:)=0.01 ; % [dB/(MHz^y cm)] 1st medium.alpha_power = 1.5; medium.BonA = 6; % create the time array t_end = 120e-6; % [s] kgrid.t_array = makeTime(kgrid, medium.sound_speed, [], t_end); % ========================================================================= % DEFINE THE SOURCE % ========================================================================= % define properties of the input signal source_strength = 10e6; % [Pa] tone_burst_freq = 2.25e6; % [Hz] tone_burst_cycles = 3; % create the input signal using toneBurst input_signal = toneBurst(1/kgrid.dt, tone_burst_freq, tone_burst_cycles); % scale the source magnitude by the source_strength divided by the % impedance (the source is assigned to the particle velocity) input_signal = (source_strength./(1000*1134)).*input_signal; % ========================================================================= % DEFINE THE ULTRASOUND TRANSDUCER % ========================================================================= % physical properties of the transducer transducer.number_elements = 16; % total number of transducer elements transducer.element_width = 1; % width of each element [grid points/voxels] transducer.element_length = 12*4/3; % length of each element [grid points/voxels] transducer.element_spacing = 0; % spacing (kerf width) between the elements [grid points/voxels] transducer.radius = inf; % radius of curvature of the transducer [m] % calculate the width of the transducer in grid points transducer_width = transducer.number_elements*transducer.element_width ... + (transducer.number_elements - 1)*transducer.element_spacing; % use this to position the transducer in the middle of the computational grid transducer.position = round([1, Ny/2 - transducer_width/2, Nz/2 - transducer.element_length/2]); % properties used to derive the beamforming delays transducer.sound_speed = 1000; % sound speed [m/s] transducer.focus_distance = 70/1000; % focus distance [m] transducer.elevation_focus_distance = 68/1000; % focus distance in the elevation plane [m] transducer.steering_angle = 0; % steering angle [degrees] % apodization transducer.transmit_apodization = 'Hanning'; transducer.receive_apodization = 'Rectangular'; % define the transducer elements that are currently active number_active_elements = 16; transducer.active_elements = ones(transducer.number_elements, 1); %transducer.active_elements = zeros(transducer.number_elements, 1); %transducer.active_elements(21:52) = 1; % append input signal used to drive the transducer transducer.input_signal = input_signal; % create the transducer using the defined settings transducer = makeTransducer(kgrid, transducer); % print out transducer properties transducer.properties; % ========================================================================= % RUN THE SIMULATION % ========================================================================= % set the input settings input_args = {'DisplayMask', transducer.active_elements_mask, ... 'PMLInside', true, 'PlotPML', true, 'PMLSize', [PML_X_SIZE, PML_Y_SIZE, PML_Z_SIZE], ... 'DataCast', DATA_CAST, 'PlotScale', [-source_strength/32, source_strength/32]}; % run the simulation %[sensor_data] = kspaceFirstOrder3D(kgrid, medium, transducer, transducer, input_args{:}); [sensor_data] = kspaceFirstOrder3D(kgrid, medium, transducer, transducer, input_args{:}); % extract a single scan line from the sensor data using the current % beamforming settings scan_line = transducer.scan_line(sensor_data); save a_scan_eun scan_line; % ========================================================================= % VISUALISATION % ========================================================================= % plot the recorded time series figure; stackedPlot(kgrid.t_array*1e6, sensor_data); xlabel('Time [\mus]'); ylabel('Transducer Element'); title('Recorded Pressure'); % plot the scan line figure; plot(kgrid.t_array*1e6, scan_line, 'k-'); xlabel('Time [\mus]'); ylabel('Pressure [au]'); title('Scan Line After Beamforming');