k-Wave User Forum » Topic: wrong Pressure of Focus in bowl focused transducer simulation

k-Wave User Forum » Topic: wrong Pressure of Focus in bowl focused transducer simulation http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation Support for the k-Wave MATLAB toolbox en-US Tue, 12 May 2026 22:02:58 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8980 Thu, 07 Dec 2023 11:26:38 +0000 Bradley Treeby 8980@http://www.k-wave.org/forum/ Are you extracting the amplitude at the fundamental and the harmonics? Note, you can avoid start-up transients in your sensor data (including the <code>p_max</code> calculation) by defining <code>sensor.record_start_index</code>. Arlen on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8978 Thu, 07 Dec 2023 06:17:34 +0000 Arlen 8978@http://www.k-wave.org/forum/ Dear Bradley Treeby, Thanks for your advice! I found that this problem happens when simulation is in nonlinear condition. And I try to use sensor.record = {'p'} and extractAmpPhase() to find the max pressure, then the results are converged. I think maybe when use sensor to record 'p_max', the sound field is unstable at the beginning, the max pressure of which is larger than when the sound field is stable, and the wrong max pressure is recorded. Bradley Treeby on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8956 Fri, 24 Nov 2023 17:33:42 +0000 Bradley Treeby 8956@http://www.k-wave.org/forum/ It sounds like the simulation hasn't converged. When you run nonlinear simulations, the grid needs to be fine enough to supported the generated higher frequency harmonics. If it's not, you will get a numerical error called spectral blocking. So in short, either reduce the amplitude of your source, or keep decreasing the grid size until the simulation output doesn't depend on the computational settings. Arlen on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8931 Mon, 09 Oct 2023 10:23:04 +0000 Arlen 8931@http://www.k-wave.org/forum/ Hi, Bradley I have changed the number of grid points (Nx, Ny, Nz), keeping the x_size (y_size, z_size) nearly same. But when I delete 'medium.BonA=7', the pressure of focus converges, changing the grid size (dx, dy, dz). Bradley Treeby on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8929 Sun, 08 Oct 2023 13:05:13 +0000 Bradley Treeby 8929@http://www.k-wave.org/forum/ Do you also change the number of grid points (Nx etc) as you change the grid spacing to keep the overall grid size constant? Arlen on "wrong Pressure of Focus in bowl focused transducer simulation" http://www.k-wave.org/forum/topic/wrong-pressure-of-focus-in-bowl-focused-transducer-simulation#post-8917 Mon, 25 Sep 2023 15:09:21 +0000 Arlen 8917@http://www.k-wave.org/forum/ Hi， I'm doing bowl focused transducer simulation, with signal of 1Mhz and 543.36e3 pa. But when I use smaller dx, dy, dz ,the pressure of focus(max pressure) is getting higher. And when I use 0.3mm of grid size, the max pressure is about 30 Mpa, and 0.2mm of grid size, the max pressure is about 40 Mpa! While the same parameters I put in focusedBowlONeil() to calculate, the max pressure at focus is just about 20 Mpa, almost the same 20 Mpa results when I use 0.5mm of grid size(dx,dy,dz). Am I use so large amplitude? My English isn't very well , hope you understand my confusion and help me to solve this problem. Thanks!!! Here are my codes: clear; clf; close all; addpath('D:\k-wave-toolbox-version-1.4\k-Wave'); Nx = 130*3; % number of grid points in the x direction Ny = 80*3; % number of grid points in the y direction Nz = 80*3; % number of grid points in the z direction dx =0.3e-3; % grid point spacing in the x direction [m] dy =0.3e-3; % grid point spacing in the y direction [m] dz =0.3e-3; % grid point spacing in the z direction [m] x_size=dx*Nx; y_size=dy*Ny; z_size=dz*Nz; sound_speed=1500; % [m/s]媒介声速 density = 1000 ; % [kg/m^3]媒介密度 % radius = 63*2; %换能器曲率半径 % diameter = 63*2 - 1; %换能器直径（奇数） kgrid = kWaveGrid(Nx, dx, Ny, dy, Nz, dz); % define the properties of the propagation medium medium.BonA=7; %水里的超声非线性参数 medium.sound_speed = sound_speed * ones(Nx, Ny, Nz); % [m/s] medium.density = density * ones(Nx, Ny, Nz); % [kg/m^3] kgrid.makeTime(medium.sound_speed,0.3,80e-6); % define a centered circular sensor % define a series of Cartesian points to collect the data sensor.mask = zeros(Nx, Ny, Nz); sensor.mask(:, :, Nz/2) = 1; % set the record mode such that only the max and peak values are stored sensor.record = {'p_max'}; % % input arguments input_args = {'PlotPML', true, ... 'DataCast', 'single', 'CartInterp', 'linear','PlotSim',true,'PlotLayout',true,'PMLInside',true,'PMLSize',[20,20,20]}; % grid_size = [Nx, Ny, Nz];%网格空间 % bowl_pos = [1, Ny/2,Nz/2];%bowl的中心点 % %radius = Nx/2;%曲率半径 % %diameter =Nx-1;%圆的直径 % focus_pos = [Nx/2, Ny/2, Nz/2 ];%与bowl_pos连线为对称轴（确定方向）和法向量有一点点像决定换能器姿态 % source.p_mask=makeBowl(grid_size, bowl_pos, radius, diameter, focus_pos, 'Plot', true); % define a time varying sinusoidal source source_mag = 543.36e3; % [Pa] % create empty array karray = kWaveArray('SinglePrecision',true); %position = [1, Ny/2, Nz/2]; position = [-x_size/2+20e-3,0,0]; % define arc properties radius = 63.2e-3; % [m] diameter = 64e-3; % [m] focus_pos = [x_size/2, 0, 0]; karray.addBowlElement(position, radius, diameter, focus_pos); source.p_mask = karray.getArrayBinaryMask(kgrid); voxelPlot(double(source.p_mask)); source_freq = 1e6; % [Hz]声源频率 sig1 = source_mag * sin(2 * pi * source_freq * kgrid.t_array); % smooth the source sig1 = filterTimeSeries(kgrid, medium, sig1,'PlotSignals',true,'PlotSpectrums',true); source_signal = zeros(1,length(sig1)); source_signal(1,1:length(sig1)) = sig1; source.p = karray.getDistributedSourceSignal(kgrid, source_signal); sensor_data = kspaceFirstOrder3DG(kgrid, medium, source, sensor, input_args{:}); sensor_data.p_max = reshape(sensor_data.p_max, [Nx, Ny]); sensor_data.p_max = reshape(sensor_data.p_max, [Nx, Ny]); figure; imagesc(kgrid.x_vec, kgrid.y_vec,sensor_data.p_max); xlabel('x [m]'); ylabel('y [m]'); %zlabel('z [m]'); title('Maximum Pressure'); c = colorbar; c.Label.String = 'pressure()'; max = sensor_data.p_max(1,1); maxi = [1,1]; [max, maxi]; for a=floor(0.3*Nx):Nx for b=1:Ny if sensor_data.p_max(a,b)>max max= sensor_data.p_max(a,b); maxi=[a,b]; end end end focus_x = 0; focusx_max = 0; focusy_max = 0; focus_y = 0; focus_l = 0; focus_r = 0; focus_u = 0; focus_d = 0; [focus_x,focus_y]; for a=floor(0.3*Nx):(Nx - 1) for b=1:(Ny - 1) if sensor_data.p_max(a,b)<0.71*max&&sensor_data.p_max(a,b+1)>0.71*max %sensor_data.p_max(a,b)=0; focus_l = b+1; end if sensor_data.p_max(a,b)>0.71*max&&sensor_data.p_max(a,b+1)<0.71*max %sensor_data.p_max(a,b)=0; focus_r = b+1; end if sensor_data.p_max(a,b)<0.71*max&&sensor_data.p_max(a+1,b)>0.71*max focus_u = a+1; end if sensor_data.p_max(a,b)>0.71*max&&sensor_data.p_max(a+1,b)<0.71*max %sensor_data.p_max(a,b)=0; focus_d = a+1; end end focus_x = focus_r - focus_l; focus_y = focus_d - focus_u; if focus_x > focusx_max focusx_max = focus_x; end if focus_y > focusy_max focusy_max = focus_y; end % if focus_d > 0 % break; % end end % define transducer parameters radius = 63.2e-3; % [m] diameter = 64e-3; % [m] frequency = 1e6; % [Hz] sound_speed = 1500; % [m/s] density = 1000; % [kg/m^3] Z = density*sound_speed; p = 543.36e3; velocity = p/Z; % [m/s] % define position vectors axial_position = 0:5e-4:150e-3; % [m] lateral_position = -15e-3:5e-4:15e-3; % [m] % evaluate pressure [p_axial, p_lateral] = focusedBowlONeil(radius, diameter, velocity, ... frequency, sound_speed, density, axial_position, lateral_position); % plot figure; subplot(4, 1, 1); plot(axial_position .* 1e3, p_axial .* 1e-6, 'k-'); xlabel('Axial Position [mm]'); ylabel('Pressure [MPa]'); subplot(4, 1, 2); plot((0:dx:(Nx-1)*dx) .* 1e3, sensor_data.p_max(:,Ny/2).* 1e-6, 'k-'); xlabel('Axial Position [mm]'); ylabel('Pressure [MPa]'); subplot(4, 1, 3); plot(lateral_position .* 1e3, p_lateral .* 1e-6, 'k-'); xlabel('Lateral Position [mm]'); ylabel('Pressure [MPa]'); subplot(4, 1, 4); plot((0:dy:(Ny-1)*dy).* sensor_data.p_max(Ny/2,:) .* 1e-6, 'k-'); xlabel('Lateral Position [mm]'); ylabel('Pressure [MPa]'); figure; imagesc(kgrid.x_vec, kgrid.y_vec,sensor_data.p_max); xlabel('x [m]'); ylabel('y [m]'); %zlabel('z [m]'); title('Maximum Pressure'); c = colorbar; c.Label.String = 'pressure()';