k-Wave User Forum » User Favorites: dsm194

k-Wave User Forum » User Favorites: dsm194 dsm194 Support for the k-Wave MATLAB toolbox en-US Wed, 13 May 2026 01:51:01 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-7232 Wed, 12 Feb 2020 16:19:48 +0000 Bradley Treeby 7232@http://www.k-wave.org/forum/ Hi DLam, I'd suggest recording a time varying pressure signal for an integer number of periods at the end of the simulation. Then spectrally decompose the signal, and sum the absorption coefficient times intensity (calculated using p^2) at each harmonic. Brad. DLam on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-7200 Tue, 28 Jan 2020 17:26:28 +0000 DLam 7200@http://www.k-wave.org/forum/ Hi Brad, I'm replying to this thread since the work here is similar to what I wanted to ask. I have been using the C++ version of the simulator, where there is no I_avg output flag available. Currently, I am using the approximation that I_avg = (P_max^2)/(2*density*soundspeed) to calculate the scalar intensity array. Can I still use Q = -div() in this case? I suppose the vector field will have to be assumed somehow, e.g. to be along the transducer beam axis. Mithun on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4849 Fri, 21 Nov 2014 16:40:12 +0000 Mithun 4849@http://www.k-wave.org/forum/ Hi Brad, That was of great help. Thank you. Kind Regards, Mithun Bradley Treeby on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4810 Sat, 01 Nov 2014 11:30:51 +0000 Bradley Treeby 4810@http://www.k-wave.org/forum/ Hi Mithun, If you plot the simulation, you'll see the wave from the edge of the aperture. It looks like a spherical wave (or cylindrical wave depending on your aperture shape), but varies in polarity on either side of planar wave front. I've created a very simple example to illustrate: <pre><code>clear all; % create the computational grid Nx = 64; % number of grid points in the x direction Ny = 32; % number of grid points in the y direction Nz = 96; % number of grid points in the z direction dx = 0.1e-3; % grid point spacing in the x direction [m] dy = 0.1e-3; % grid point spacing in the y direction [m] dz = 0.1e-3; % grid point spacing in the z direction [m] kgrid = makeGrid(Nx, dx, Ny, dy, Nz, dz); % define the properties of the propagation medium medium.sound_speed = 1500; % create time array kgrid.t_array = makeTime(kgrid, medium.sound_speed, 0.5, 4e-6); % create initial pressure distribution source.p0 = zeros(Nx, Ny, Nz); source.p0(20, 14:19, 38:59) = 1; % define sensor mask sensor.mask = zeros(Nx, Ny, Nz); sensor.mask(Nx/2, Ny/2, Nz/4) = 1; sensor.mask(3*Nx/4, Ny/2, Nz/2) = 1; % input arguments input_args = {'PlotPML', false, 'DataCast', 'single'}; % run the simulation sensor_data = kspaceFirstOrder3D(kgrid, medium, source, sensor, input_args{:}); % plot the simulated sensor data figure; plot(sensor_data.') ylabel('Pressure'); xlabel('Time Step'); legend('Off Axis', 'On Axis');</code></pre> Hope that helps, Brad. Mithun on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4789 Fri, 24 Oct 2014 10:33:07 +0000 Mithun 4789@http://www.k-wave.org/forum/ Hi Brad, Thank you for the valuable comments. Yes of-course, I can zero the first part of the signal. But I would like to understand whats happening actually. The high signal is coming not in the very beginning, but at around 100th sample for the middle element. It is a bipolar signal with negative cycle first and positive next. Then I thought it may be a reflection with phase shift of input signal with 1 cycle. How can I recognize an edge wave? Combining sensor data is really useful. Thank you. Regards, Mithun Bradley Treeby on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4784 Tue, 21 Oct 2014 16:04:57 +0000 Bradley Treeby 4784@http://www.k-wave.org/forum/ Hi Mithun, This function was added as part of k-Wave V1.1 (released a couple of weeks ago). The signal at the beginning will be your input signal. Inside the code, the source and sensor act independently, so if they are in the same position, the sensor will also record the input signal. There's a couple of ways you could remove this: (1) You could zero the first part of the signal, e.g., following the "Simulating B-mode Ultrasound Images Example". (2) You could run a separate simulation in a homogeneous medium, and subtract the recorded signal from the signal you record when you have scatterers. Are you sure it's a reflection, and not an edge wave from the edge of the transducer? Brad. Mithun on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4783 Tue, 21 Oct 2014 15:02:38 +0000 Mithun 4783@http://www.k-wave.org/forum/ Dear Brad, I fixed it by adding the code you provided in another post to KwaveTransducer class. Thats amazing, it saved lot of time. Thank you. My problem is that I am getting a huge initial signal in all elements. This is a phase inverted bipolar signal. Is it the input signal itself? How can I get rid of it? My transducer is below PML I guess. Even after this initial signal, I have another weak reflection also from unexpected depth. Any ideas on it looking at the code above? Thank You, Kind Regards, Mithun Mithun on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4781 Tue, 21 Oct 2014 14:08:29 +0000 Mithun 4781@http://www.k-wave.org/forum/ Dear Brad, Thanks for the response. Great to know about apodization in c++ code. I am not able to use transducer.combine_sensor_data(sensor_data). Its throwing an error "No appropriate method, property, or field combine_sensor_data for class kWaveTransducer". Do I have to write the method explicitly? Please help me on this. I am not going to use scan_line method as I am interested in data from each elements separately. Thank You, Regards, Mithun Bradley Treeby on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4771 Thu, 16 Oct 2014 11:40:56 +0000 Bradley Treeby 4771@http://www.k-wave.org/forum/ Hi Mithun, There is currently a difference in the way the C++ and MATLAB codes return the data when using an object of the <code>kWaveTransducer</code> class as a sensor. The MATLAB code returns one time series per active transducer element, while the C++ code returns one time series per grid point that forms part of the active transducer. If using the C++ code, before using the <code>scan_line</code> method, you should call: <code>sensor_data = transducer.combine_sensor_data(sensor_data);</code> This will combine the sensor data to give a single time series per active transducer element, analogous to the output from the MATLAB code. Note that transmit apodization is not currently implemented in the C++ code, i.e., the apodization is always the same as: <code>transducer.transmit_apodization = 'Rectangular';</code> (see the discussion <a href="http://www.k-wave.org/forum/topic/question-about-the-transducercombine_sensor_data-method">here</a>). Let me know if that doesn't solve the problem. Brad. Mithun on "Image Artefacts - Ultrasound Imaging" http://www.k-wave.org/forum/topic/image-artefacts-ultrasound-imaging#post-4767 Sun, 12 Oct 2014 09:34:55 +0000 Mithun 4767@http://www.k-wave.org/forum/ Dear All, I am trying to do a simple US image simulation using K-wave [Transducer as source and sensor]. Quite new to it. I created a medium and then a ball in the middle with different speed of sound and density. I expected only the reflection from the ball. Instead I am seeing a very high signal in the beginning of time trace (May be input signal?, How to get rid?) and different reverberations at almost all depths. What may be the reason for this? I put my transducer below PML. I am pasting my script below. Can anyone please throw your views on it? It would be of great help for me to move forward. clear all; clc % simulation settings DATA_CAST = 'single'; % ========================================================================= % DEFINE THE K-WAVE GRID % ========================================================================= % set the size of the perfectly matched layer (PML) PML_X_SIZE = 20; % [grid points] PML_Y_SIZE = 10; % [grid points] PML_Z_SIZE = 10; % [grid points] % set total number of grid points not including the PML Nx = 512 - 2*PML_X_SIZE; % [grid points] Ny = 512 - 2*PML_Y_SIZE; % [grid points] Nz = 128 - 2*PML_Z_SIZE; % [grid points] % set desired grid size in the x-direction not including the PML x = 30e-3; % [m] % calculate the spacing between the grid points dx = x/Nx; % [m] dy = dx; % [m] dz = dx; % [m] % create the k-space grid kgrid = makeGrid(Nx, dx, Ny, dy, Nz, dz); % Total time, 2 times for US t_end = (x/1500)*2; % [s] [kgrid.t_array,dt] = makeTime(kgrid, 1500, [], t_end); % ========================================================================= % DEFINE THE INPUT SIGNAL % ========================================================================= % define properties of the input signal source_strength = 5e6; % [Pa] tone_burst_freq = 1.5e6; % [Hz] tone_burst_cycles = 2; % 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./(1500*1000)).*input_signal; medium.sound_speed =1500; % [m/s] medium.density = 1000; % [kg/m^3] % medium.alpha_coeff = 0.0022; % [dB/(MHz^y cm)] % medium.alpha_power = 2; % medium.BonA = 6; % define properties density_map = 1000*ones(Nx_tot, Ny_tot, Nz_tot); medium.density = density_map (:,:,:); ball_radius = round(.1e-3/dx); % [mm] p1 = 962*makeBall(Nx, Ny, Nz,PML_X_SIZE + round(16e-3/dx),PML_Y_SIZE + round(15e-3/dx), PML_Z_SIZE + Nz/2, ball_radius); medium.density = medium.density + p1; % DEFINE THE ULTRASOUND TRANSDUCER % ========================================================================= % physical properties of the transducer transducer.number_elements = 64; % total number of transducer elements transducer.element_width = 2; % width of each element [grid points/voxels] transducer.element_length = 24; % 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([PML_X_SIZE+1, Ny/2 - transducer_width/2, Nz/2 - transducer.element_length/2]); % properties used to derive the beamforming delays transducer.sound_speed = 1500; % sound speed [m/s] transducer.focus_distance = 15e-3; % focus distance [m] transducer.elevation_focus_distance = 16e-3 + PML_X_SIZE; % 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 transducer.active_elements = zeros(transducer.number_elements, 1); transducer.active_elements(1:transducer.number_elements) = 1; 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 % ========================================================================= sensor.mask=transducer.all_elements_mask; % set the input settings input_args = {... 'PMLInside', false, 'PlotSim', false, 'PMLSize', [PML_X_SIZE, PML_Y_SIZE, PML_Z_SIZE], ... 'DataCast', DATA_CAST}; % run the simulation sensor_data = kspaceFirstOrder3DC(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); % reshape sensor data to y, z, t sensor_data_rs = reshape(sensor_data, transducer.number_elements*transducer.element_width, transducer.element_length, length(kgrid.t_array)); Thank You Kind Regards, Mithun Bradley Treeby on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-4096 Wed, 20 Nov 2013 10:24:38 +0000 Bradley Treeby 4096@http://www.k-wave.org/forum/ Hi Anthony, Thanks for your code - I see what you mean. You could try using the <code>gradientSpect</code> function in k-Wave to see if that makes a difference to your calculation of the divergence. This computes the derivative in the same way as within the simulation functions. In regards to the accuracy of the intensity calculation, I've added some code to shift the velocity to the non-staggered spatial and temporal grid before multiplying by the pressure, and this seems to improve things. This will go into the V1.1 release. Brad. Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-4092 Fri, 15 Nov 2013 15:45:57 +0000 Anthony 4092@http://www.k-wave.org/forum/ This small piece of code illustrates what I have just written (you just have to specify the correct path in the 3 first lines). Anthony ********************************* CODE ****************************** clear all close all clc infilename = 'C:\Users\Grisey\Desktop\checkSource\input.h5'; outfilename = 'C:\Users\Grisey\Desktop\checkSource\output.h5'; pathToKwave = '"C:\Users\Grisey\k-wave\k-wave-toolbox-version-1.0-cpp-windows-binaries\Windows Binaries\kspaceFirstOrder3D-OMP.exe"'; %(don't forget the ") %% Paramètres myPMLsize = 20; medium.density = 1000; medium.sound_speed = 1500; source_mag = 1e5; %% Grille de calcul f = 1e5; N = 128; h = 1e-3; kgrid = makeGrid(N, h, N, h, N, h); %% Definition du transducteur R = 10e-3; transducerMask = zeros(N,N,N); xgrid = h*((1:N)-N/2); ygrid = h*((1:N)-N/2); zgrid = h*((1:N)-N/2); [X,Y,Z] = meshgrid(xgrid,ygrid,zgrid); Rmap = (X.^2+Y.^2+Z.^2).^0.5; transducerMask(Rmap<=R)=1; %% Simulation % create the time array [kgrid.t_array, dt] = makeTime(kgrid, medium.sound_speed); periode = 1/f; recStartTime = numel(kgrid.t_array)-2*round(periode/dt)+1; source.p_mask = logical(transducerMask); source.p = source_mag*sin(2*pi*f*kgrid.t_array); source.p_mode = 'dirichlet'; % assign the input options input_args = {'PlotPML', true, 'PMLSize', myPMLsize,'Smooth',[false true true]}; % run the simulation toto = ones(N, N, N); sensor.mask = logical(toto); clear toto; kspaceFirstOrder3D(kgrid, medium, source, sensor, input_args{:},'SaveToDisk',infilename); status = system([pathToKwave ' -i ' infilename ' -o ' outfilename ' --p_rms --I_avg -s' num2str(recStartTime)]); %% Lecture du résultat p_rms = reshape(h5read(outfilename,'/p_rms'),N,N,N); Ix_avg = reshape(h5read(outfilename,'/Ix_avg'),N,N,N); Iy_avg = reshape(h5read(outfilename,'/Iy_avg'),N,N,N); Iz_avg = reshape(h5read(outfilename,'/Iz_avg'),N,N,N); H = -divergence(Ix_avg,Iy_avg,Iz_avg)/h; % Calcul de Ir Ianalytical = source_mag^2*R^2./(2*medium.sound_speed*medium.density.*Rmap.^2); theta = acos(Z./Rmap); phi = atan(Y./X)+pi*(X<0); Ix_analytical = Ianalytical.*cos(phi).*sin(theta); Iy_analytical = Ianalytical.*sin(phi).*sin(theta); Iz_analytical = Ianalytical.*cos(theta); figure(1) subplot(1,3,1) imagesc(squeeze(H(round(N/2),:,:)).*(1-squeeze(transducerMask(round(N/2),:,:)))) caxis([-0.1e5 2e4]) title('-div(k-wave intensity) (order 1 uncentered divergence)') subplot(1,3,2) Hanalytical = -divergence(Ix_analytical,Iy_analytical,Iz_analytical)/h; imagesc(squeeze(Hanalytical(round(N/2),:,:)).*(1-squeeze(transducerMask(round(N/2),:,:)))) caxis([-0.1e5 2e4]) title('-div(analytical intensity) (order 1 uncentered divergence)') subplot(1,3,3) divx = (Ix_avg(3:end,2:end-1,2:end-1)-Ix_avg(1:end-2,2:end-1,2:end-1))/2/h; divy = (Iy_avg(2:end-1,3:end,2:end-1)-Iy_avg(2:end-1,1:end-2,2:end-1))/2/h; divz = (Iz_avg(2:end-1,2:end-1,3:end)-Iz_avg(2:end-1,2:end-1,1:end-2))/2/h; HnumericalOrder2 = -divx-divy-divz; imagesc(squeeze(HnumericalOrder2 (round(N/2),:,:)).*(1-squeeze(transducerMask(round(N/2),2:end-1,2:end-1)))) caxis([-0.1e5 2e4]) title('-div(k-wave intensity) (order 2 centered divergence)') Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-4091 Fri, 15 Nov 2013 13:54:18 +0000 Anthony 4091@http://www.k-wave.org/forum/ Hi Brad, According to some tests I have made, I think it comes from the numerical calculation of the divergence. In the simple case of a radiating sphere, with a linear medium and no attenuation, -div(pu)=0 theoretically. If I compute it numerically with an uncentered scheme of order 1, I get strong positive and negative zones along the axis of the cartesian grid. If I use a 2 order centered scheme, it becomes much better. However, it does not seem sufficient to solve the problem in the case of a focused field... Is it possible for you to try to do the calculation of the divergence in the spectral domain in k-wave in order to see if it solves the problem? Best regards, Anthony Bradley Treeby on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-4047 Thu, 31 Oct 2013 20:12:43 +0000 Bradley Treeby 4047@http://www.k-wave.org/forum/ Hi Anthony, We are still looking into the heat deposition issue. We've had a few blackboard sessions, but I think more questions have been raised than answered! Certainly you wouldn't expect to see a heat sink, but exactly how the heat deposition should be computed depends partly on where the energy loss terms appear in the conservation equations. If we get any closer, I'll let you know. Regarding the variation in the expected power, this might be due to approximations in the way intensity is calculated in V1.0 (see <a href="http://www.k-wave.org/forum/topic/defining-incident-acoustic-plane-wave-intensity#post-3723">here</a>). We've changed the way this is calculated ready for V1.1. Brad. Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-4044 Thu, 31 Oct 2013 17:35:17 +0000 Anthony 4044@http://www.k-wave.org/forum/ Hi Brad, I must admit that I couldn't solve the mystery (still using 2*alpha*I in the meantime). Did you have some new insight about this? I wrote a piece of code that shows the problem, if someone has an idea... By the way, the acoustic power I get by integrating intensity over a plane is smaller than the value I should theoretically get : I guess it comes from the way I compute the amplitude of the velocity source : u0 = sqrt( 2 * acousticPower / ( rho * soundSpeed * transducerSurface ) ) Is there another way to compute it (without the plane wave approximation, maybe...)? Best regards, Anthony ******************************** SAMPLE CODE ************************************ clear all close all clc %% Files work_dir = uigetdir('Set working directory (where the files will be saved) '); kspaceFirstOrder3D_OMP_path = uigetdir('Where is kspaceFirstOrder3D-OMP.exe?'); acousticSimulationInputFileName = 'input.h5'; acousticSimulationOutputFileName = 'output.h5'; %% Parameters % Source acousticPower = 80; % [W] source_freq = 2.5e5; %[Hz] % Physical parameters rho = 1000; soundSpeed = 1500; %% Geometry R = 40e-3; % [m] Rap = 30e-3; % [m] h = 0.5e-3; % [m] N = ceil(62e-3/h)+20; % acoustic simulation grid size %% k-wave parameters % Domain geometry myPMLsize = 10; % Source transducerSurface = 2*pi*R*(R-sqrt(R^2-Rap^2)); source_mag = sqrt(2 * acousticPower / (rho * soundSpeed * transducerSurface)); %velocity source kgrid = makeGrid(N, h, N, h, N, h); % define a curved transducer element xCenter = round(N/2); yCenter = round(N/2); zCenter = round(R/h)+2; transducerMask = zeros(N,N,N); xGrid = (1:N)'*ones(1,N); yGrid = ones(N,1)*(1:N); for zmap=1:(zCenter-(sqrt(R^2-Rap^2)/h)); rmap2D = round(((xGrid-xCenter).^2+(yGrid-yCenter).^2+(zmap-zCenter).^2).^0.5); transducerMask(:,:,zmap) = 0+(rmap2D==round((R)/h)); end clear xGrid yGrid rmap2D % voxelPlot(transducerMask) % define the properties of the propagation medium medium.density = 1000; % [kg/m3] medium.sound_speed = 1500; % [m/s] % create the time array [kgrid.t_array, dt] = makeTime(kgrid, medium.sound_speed); periode = 1/source_freq; recStartTime = numel(kgrid.t_array)-round(periode/dt)+1; source.u_mask = logical(transducerMask); % define a time varying sinusoidal source [I,J,K] = ind2sub([N,N,N],find(source.u_mask)); I = xCenter - I; J = yCenter - J; K = zCenter - K; distCenter = (I.^2+J.^2+K.^2).^0.5; timeSine = sin(2*pi*source_freq*kgrid.t_array); source_ux0 = abs(I./distCenter)*source_mag; source_uy0 = abs(J./distCenter)*source_mag; source_uz0 = abs(K./distCenter)*source_mag; source.ux = source_ux0*timeSine; source.uy = source_uy0*timeSine; source.uz = source_uz0*timeSine; % assign the input options input_args = {'PlotPML', true, 'PMLSize', myPMLsize,'Smooth',[false true true]}; % run the simulation sensor.mask = logical(ones(N, N, N)); kspaceFirstOrder3D(kgrid, medium, source, sensor, input_args{:},'SaveToDisk',[work_dir acousticSimulationInputFileName]); clear kgrid medium source sensor status = system(['"' kspaceFirstOrder3D_OMP_path '\kspaceFirstOrder3D-OMP.exe" -i ' [work_dir acousticSimulationInputFileName] ' -o ' [work_dir acousticSimulationOutputFileName] ' --I_avg -s' num2str(recStartTime)]); %% Load the results Ix_avg = h5read([work_dir acousticSimulationOutputFileName],'/Ix_avg'); Ix_avg = reshape(Ix_avg,N,N,N); Iy_avg = h5read([work_dir acousticSimulationOutputFileName],'/Iy_avg'); Iy_avg = reshape(Iy_avg,N,N,N); Iz_avg = h5read([work_dir acousticSimulationOutputFileName],'/Iz_avg'); Iz_avg = reshape(Iz_avg,N,N,N); %% Compute heat source term heatSource = -divergence(Ix_avg,Iy_avg,Iz_avg)/h; figure(1) imagesc(squeeze(heatSource(:,round(N/2),:))') xlabel('X') ylabel('Z') title('Heat source term : -div(I_{avg})') % Acoustic power Pac_from_I = max(squeeze(sum(sum(Iz_avg(myPMLsize:end-myPMLsize,myPMLsize:end-myPMLsize,myPMLsize:end-myPMLsize),2),1))*h^2); disp(['Acoustic power : ' num2str(Pac_from_I) 'W']) Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-3829 Fri, 26 Jul 2013 09:27:48 +0000 Anthony 3829@http://www.k-wave.org/forum/ Bradley Treeby on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-3825 Thu, 25 Jul 2013 18:46:45 +0000 Bradley Treeby 3825@http://www.k-wave.org/forum/ Hi Anthony, Interesting question! Leave it with us for a short time, I'd like to check a few things before I give you an answer. Brad. Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-3818 Thu, 25 Jul 2013 14:19:15 +0000 Anthony 3818@http://www.k-wave.org/forum/ I realize my question is not clear. I understand that the formula simply comes from the Green formula, but... when I apply it to a focused field, I get a strong heat source right before the focal point and a strong (but weaker due to absorption) heat sink right after the focal point, what results in a strange thermal behaviour. (<a href="http://www.casimages.com/img.php?i=130725031106339995.jpg" title="example">see example</a>) From your experience, is this normal? I must probably miss some point... Regards, Anthony on "Using intensity as heat source term (nonlinear acoustics)" http://www.k-wave.org/forum/topic/using-intensity-as-heat-source-term-nonlinear-acoustics#post-3802 Tue, 23 Jul 2013 08:21:37 +0000 Anthony 3802@http://www.k-wave.org/forum/ Hi Bradley, Just to be sure, when you write Q = - div(), the source term become negative at some places, is it normal? If yes, how do you physically explain it? Anthony PS : I am sorry, I suppose it is probably in Pierce's book but I don't have it and there is nothing in the book of Hamilton & Blackstock... Do you know another ref about this particular topic? jrosenfield on "Defining Incident Acoustic Plane Wave Intensity" http://www.k-wave.org/forum/topic/defining-incident-acoustic-plane-wave-intensity#post-3726 Tue, 18 Jun 2013 04:21:25 +0000 jrosenfield 3726@http://www.k-wave.org/forum/ Hi Brad, I also meant to let you know that when I compute the average squared acoustic pressure using the time series of pressure values recorded by k-Wave, the average that I calculate also seems susceptible to the same frequency dependence as the average acoustic intensity. Specifically, the average squared pressure approaches the theoretical average as I decrease the source frequency away from the maximum frequency supported by the grid. I am not sure why this would be the case, since my understanding is that the pressure itself is not recorded on a staggered grid, so the issue causing problems with the average intensity calculation should not be a problem with pressure as you mentioned earlier. -Jon