k-Wave User Forum » User Favorites: Pezeril

k-Wave User Forum » User Favorites: Pezeril pezeril Support for the k-Wave MATLAB toolbox en-US Wed, 13 May 2026 01:04:40 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-7007 Thu, 22 Aug 2019 15:44:23 +0000 Bradley Treeby 7007@http://www.k-wave.org/forum/ Hi Patrick, Unfortunately, it is still not in the current release. Yes, that looks correct. Brad. Patrick on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6995 Sun, 04 Aug 2019 05:25:55 +0000 Patrick 6995@http://www.k-wave.org/forum/ As a follow-up. I modified the 3D elastic code to record the normal components of shear stress. Am I correct in using <code>sxx = sxx_split_x + sxx_split_y + sxx_split_z</code>? (and similarly for syy and szz and of course with the correct indexing) Thanks again, PC Patrick on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6994 Sun, 04 Aug 2019 04:25:24 +0000 Patrick 6994@http://www.k-wave.org/forum/ Hello, Has the ability to record the components of shear stress been implemented in the current release? If so, how can this be done? Thanks, PC svch1mail on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6802 Wed, 27 Mar 2019 15:48:28 +0000 svch1mail 6802@http://www.k-wave.org/forum/ Hi, I had the same problem. Thank you Brad for explanations, they work, but additionally in pstdElastic2D.m it need to add after 1017 line someting like: SS.p=sensor_data.p; SS.sxy=sensor_data.sxy; sensor_data = SS;%sensor_data.p; Thanks, Sergiy ad000000 on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6401 Mon, 09 Apr 2018 19:34:55 +0000 ad000000 6401@http://www.k-wave.org/forum/ Thanks Brad! Bradley Treeby on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6396 Sun, 08 Apr 2018 18:44:39 +0000 Bradley Treeby 6396@http://www.k-wave.org/forum/ Hi Andy, For k-Wave V1.2, to hack it into the <code>pstdElastic2D</code> in the simplest way, immediately before the time loop (line 617) add: <code>sensor_data.sxy = castZeros([num_sensor_points, num_recorded_time_points]);</code> Then after line 814, add the following: <code>sensor_data.sxy(:, file_index) = sxy_split_x(sensor_mask_index) + sxy_split_y(sensor_mask_index);</code> Similarly for 3D. Note, I haven't tried this, so let me know if it doesn't work for some reason. Brad. Bradley Treeby on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6393 Sun, 08 Apr 2018 18:27:36 +0000 Bradley Treeby 6393@http://www.k-wave.org/forum/ Hi ad000000, The dispersion is due to the low-order numerical schemes typically used in FEM. Switching to COMSOL is not likely to solve this issue. For a scatterer in an otherwise homogeneous domain, a combined FEM-BEM approach could work well. There are lots of open source FEM projects (e.g., <a href="https://fenicsproject.org">FEniCS</a>), but I haven't really used them enough to make a recommendation. Brad. ad000000 on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6374 Mon, 26 Mar 2018 18:24:04 +0000 ad000000 6374@http://www.k-wave.org/forum/ Hi Brad, Could you give a little bit more advice on the above hack? I'm modeling the typical NDT case of a longitudinal wave launched in an acrylic wedge, which gets converted to shear wave in a steel part. I'd like to extract the pressure (or velocity) from the shear wave in the steel, but there is also a longitudinal wave there (a couple, actually). So the "u" or "u_max" data would map both waves together. The shear stress, if it could be extracted and saved, would solve that problem. Thanks, Andy ad000000 on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6373 Mon, 26 Mar 2018 18:14:23 +0000 ad000000 6373@http://www.k-wave.org/forum/ Thanks Brad. I really like K-wave, for the versatility. I'll try the convergence test. However I take your point about FEM being a better choice for elastic circles, large discontinuities in material properties etc. Can you recommend a FEM program for this? I know that PZflex has a core solver that is susceptible to considerable numerical dispersion over distances on the order of hundreds of wavelengths. I don't know if Comsol has same problem, but it's a bit pricey... Bradley Treeby on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6363 Wed, 21 Mar 2018 21:43:21 +0000 Bradley Treeby 6363@http://www.k-wave.org/forum/ Hi ad000000, Glad you managed to solve most of the issues. In general, modelling a circle in an elastic medium with a very large impedance mismatch is a challenging problem for the type of numerical model used in k-Wave (FEM would handle this much better). There are two issues. The first is staircasing, which effects elastic simulations more than fluid simulations (see <a href="http://bug.medphys.ucl.ac.uk/papers/2017-Robertson-JASA.pdf">this paper</a> for some details). The second is the large discontinuity in medium properties. Both of these can be solved by increasing the sampling. If you want to gain some confidence in your numerical results, I would suggest a convergence test, i.e., reduce dx and dt, and check you get the same answer. Brad. ad000000 on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6359 Wed, 21 Mar 2018 02:18:32 +0000 ad000000 6359@http://www.k-wave.org/forum/ ...continued: 1. Assigning uy=signal, ux=0 is NOT the way to make a shear wave. Need to use sxy. 2. I was able to set the hole material to very low density and pretty low velocity. Max attenuation I could set was pretty low, at about 20 dB/cm. This still gave substantial oscillations withing the hole, but I don't *think* they altered the scattering of the main wave much. Thus I think I got a reasonable approximation for scattering from a hole. I also tried using Dirichlet BCs and setting the hole as a source with u=0, but this only worked for a longitudinal incident wave, since I could not use the "u" source to make a shear wave, and could not use Dirichlet BC for "s" sources. Setting u=0 actually models an infinitely hard "hole", rather than air. Even so, the program still crashed. Apparently using the Dirichlet approach to make hard boundaries requires much care. ad000000 on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6358 Wed, 21 Mar 2018 02:04:15 +0000 ad000000 6358@http://www.k-wave.org/forum/ Follow-up: I solved these issues, mostly. ad000000 on "model shear wave scattering from hole" http://www.k-wave.org/forum/topic/model-shear-wave-scattering-from-hole#post-6354 Tue, 20 Mar 2018 05:23:45 +0000 ad000000 6354@http://www.k-wave.org/forum/ Hello, I'm trying to model scattering from a plane shear wave (signal=wavelet) from a hole. My material is steel, and I set the hole material velocity and density to be less than steel. I get several problems: 1. I can make a compression plane wave emanating from the left edge of the grid, by setting source.ux=0 and source.uy=signal. But if I reverse them, to make a shear wave, it gives nonsense results. 2. To simulate a hole (air inside), I set the density very low. However the "hole" seems to generate its own energy. I tried to kill that by setting attenuation high, but the program ground to a halt. Thanks for any thoughts! Code: clear source sensor medium kgrid % define variables freq=5e6; epw=12; cp=5850; cs=3230; rho=8030; SDHdia=1.5e-3; %==define simulation area (remember, horizontal is Y direction)================== xlen=20e-3; ylen=16e-3; % ==================== setup grid ========================================= lambda=cp/freq; dx=lambda/epw; dy=lambda/epw; Nx= ceil(xlen/dx); Ny= ceil(ylen/dy); kgrid = makeGrid(Nx, dx, Ny, dy); %====set up material regions========================================= medium.sound_speed_compression =ones(kgrid.Nx, kgrid.Ny)*cp; medium.sound_speed_shear =ones(kgrid.Nx, kgrid.Ny)*cs; medium.alpha_coeff_compression = ones(kgrid.Nx, kgrid.Ny)*.03/(freq/1e6); %to give dB/(MHz^2 cm) medium.alpha_coeff_shear = ones(kgrid.Nx, kgrid.Ny)*.03/(freq/1e6)*4; medium.density=ones(kgrid.Nx, kgrid.Ny)*rho; SDHloc(1)=0; SDHloc(2)=2e-3; SDHdia=SDHdia/dx; [SDHcenterGridx,SDHcenterGridy]=cartPt2Grid(kgrid,[SDHloc(1) SDHloc(2)]); disc = makeDisc(Nx, Ny, SDHcenterGridx, SDHcenterGridy, SDHdia/2); medium.sound_speed_compression(disc == 1) = cp*.9; medium.sound_speed_shear(disc == 1) = cs*.9; medium.density(disc==1)=60; medium.alpha_coeff_compression(disc==1) = 1/(freq/1e6); %to give dB/(MHz^2 cm) medium.alpha_coeff_shear(disc==1) = 4/(freq/1e6)*4; t_end=.3*2*ylen/cs; %if ReverseFire==1; t_end=t_end/2;end dt=.3*dx/(cp); %set time steps (C should be Cmax) kgrid.t_array = 0:dt:t_end; Nt=length(kgrid.t_array); PML=20*dx; %======================= setup source =================================================== signal = wavelet(freq,4,0,2/freq,kgrid.t_array); source.u_mask = zeros(kgrid.Nx,kgrid.Ny); % or source.p_mask source.u_mask(:,1)=1; source.ux=signal/(rho*cp); source.uy=0; %====setup sensor ===================================================== sensor.mask=zeros(kgrid.Nx, kgrid.Ny); %sensor.mask(1:5:kgrid.Nx,1:5:kgrid.Ny)=1; sensor.mask(1:2:kgrid.Nx,1:2:kgrid.Ny)=1; sensor.record = {'u_max_all','u'}; %{'p','p_max'}; %===== display_mask ========================================= display_mask=disc; %==================== RUN ================================================= figpos=get(0,'DefaultFigurePosition'); set(0,'DefaultFigurePosition',[340,275,980,735]); input_args = {'DisplayMask', display_mask, 'PMLInside', true, 'PlotPML', false,'PlotScale','auto','DataCast', 'single'};%,'RecordMovie',true}; sensor_data = pstdElastic2D(kgrid, medium, source, sensor, input_args{:}); %sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:}); set(0,'DefaultFigurePosition',figpos); kmodel=sensor_data; clear sensor_data; kmodel.freq=freq; kmodel.Nx=Nx;kmodel.Ny=Ny; kmodel.x_vec=kgrid.x_vec;kmodel.y_vec=kgrid.y_vec; kmodel.Nt=length(kgrid.t_array); kmodel.epw=epw; save kmodel.mat kmodel -v7.3 uall=sqrt(kmodel.ux_max_all.^2+kmodel.uy_max_all.^2); u=sign(kmodel.ux.*kmodel.uy).*sign(kmodel.ux).*sqrt(kmodel.ux.^2+kmodel.uy.^2); %uh=abs(hilbert(u'));uh=uh'; %nox=ceil(kmodel.Nx/5);noy=ceil(kmodel.Ny/5); nox=ceil(kmodel.Nx/2);noy=ceil(kmodel.Ny/2); dat=zeros(nox,noy);%may need to reverse this is Nx, Ny are in wrong order for t=1:size(u,2); for j=1:nox for k=1:noy; %dat(j,k,t)=uh((k-1)*nox+j,t); dat(j,k,t)=u((k-1)*nox+j,t); %daty(j,k,t)=kmodel.uy((k-1)*nox+j,t); end end end % Record the movie %maxdat=max(max(max(dat))); %ax=[0 maxdat/500]; nframes=50; step=floor((kmodel.Nt-1)/nframes); figure;imagesc(kmodel.y_vec,kmodel.x_vec,(dat(:,:,200)));axis equal;%caxis(ax) set(gca,'nextplot','replacechildren'); ax=caxis; for j= 1:nframes; t=1+(j-1)*step; imagesc(kmodel.y_vec,kmodel.x_vec,(dat(:,:,t)));caxis(ax) F(j) = getframe; end % Play the movie movie(F,1,1) Bradley Treeby on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6212 Fri, 10 Nov 2017 15:22:49 +0000 Bradley Treeby 6212@http://www.k-wave.org/forum/ Hi Lenny, You can define which parameters you would like to record using the <code>sensor.record</code> input (see the documentation for details). Unfortunately, it's not currently possible to record the components of the stress (only the particle velocity or pressure). If you're comfortable with MATLAB, it is pretty easy to hack into the code (add it after line 814 in V1.2). It's on the list for a future release. Brad. xzy9588 on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6202 Wed, 08 Nov 2017 19:34:43 +0000 xzy9588 6202@http://www.k-wave.org/forum/ This is my code clear all close all % create the computational grid Nx = 100; % number of grid points in the x (row) direction Ny = 200; % number of grid points in the y (column) 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] kgrid = kWaveGrid(Nx, dx, Ny, dy); % define the properties of the upper layer of the propagation medium medium.sound_speed_compression = 1500 * ones(Nx, Ny); % [m/s] medium.sound_speed_shear = zeros(Nx, Ny); % [m/s] medium.density = 1000 * ones(Nx, Ny); % [kg/m^3] % define the properties of the lower layer of the propagation medium medium.sound_speed_compression(Nx/2:end, :) = 2000; % [m/s] medium.sound_speed_shear(Nx/2:end, :) = 800; % [m/s] medium.density(Nx/2:end, :) = 1200; % [kg/m^3] % create the time array cfl = 0.05; % Courant-Friedrichs-Lewy number t_end = 5e-6; % [s] kgrid.makeTime(max(medium.sound_speed_compression(:)), cfl, t_end); % define a source source.s_mask = makeLine(Nx, Ny, [Nx/2,1],[1,Ny/2]); %source.u_mask = makeLine(Nx, Ny, [Nx/2,1],[1,Ny/2]); %create initial pressure distribution using makeDisc disc_magnitude = 400; % [Pa] source.p0 = disc_magnitude*source.s_mask ; % define a time varying source %source_freq = 1e6; % [Hz] %source_mag = 0.00002; %source.ux = source_mag * sin(2 * pi * source_freq * kgrid.t_array); %source.uy = source_mag * sin(2 * pi * source_freq * kgrid.t_array); %display_mask=source.u_mask; sensor.mask=makeLine(Nx,Ny,[40,40],[60,60]); % define input arguments input_args = {'PlotScale', [-100, 100, -70, 70], 'PlotPML', false,... 'DisplayMask', sensor.mask, 'DataCast', 'single'}; % run the simulation sensor_data = pstdElastic2D(kgrid, medium, source, sensor, input_args{:}); imagesc(sensor_data); colormap(getColorMap); colorbar; xzy9588 on "Failed to record shear stress." http://www.k-wave.org/forum/topic/failed-to-record-shear-stress#post-6201 Wed, 08 Nov 2017 18:56:19 +0000 xzy9588 6201@http://www.k-wave.org/forum/ Hello everyone, I tried to record the magnitude of shear stress using linear sensor mask to see how longitudinal wave can generate shear wave in a layered medium. But we failed to record the shear stress in sensors. How can I record the magnitude of shear stress as well as the pressure of longitudinal wave? Thank you Lenny