Hi, I have been using this toolbox for a few months and have had little problems although I am now a bit stuck. I have worked through the toolbox manual and I am still unsure. As I have had it working in the past, I feel like I am missing something obvious but I just can't figure it out.
I am using the toolbox to simulate an ultrasonic transducer measuring the film thickness between a piston ring and cylinder liner, for example see the linked paper. (Piezoelectric sensors to monitor lubricant film thickness at piston–cylinder contacts in a fired engine, Mills 2012 https://www.researchgate.net/publication/258178367_Piezoelectric_sensors_to_monitor_lubricant_film_thickness_at_piston-cylinder_contacts_in_a_fired_engine).
So this involves running the simulation a series of times whilst moving the "piston ring" across the domain. Within each simulation the transducer pulses a 5MHz wave, which then reflects from the far side of the liner and is recorded by the same transducer, see Figure 5 in the refered paper. The simulation runs ~1000 times for the piston ring to move across the domain, which is quite computationally expensive but I have access to a supercomputer and so that aspect isn't a problem.
Previous simulations I have run have worked fine and produced results like those I expected, (a decrease in reflected signal for the piston ring being aligned with the transducer) although I noticed an error in my input signal, so the actual input was close to 100MHz instead of 5MHz. After correcting this, I have been getting incorrect trends in results, e.g. when the piston ring is over the transducer I am getting a greater reflection when I should be getting less. The only part of the code that has been changed is the input and I am struggling to figure out why this has led to this.
Due to this I have a few questions that I am hoping you can help me with:
1. Is there any reason why this toolbox should not be able to simulate systems that have thin oil films in? (~up to 30um).
2. Is there any reason why changing the input frequency would cause might output to produce the opposite trend that it was previously producing?
3. From looking at my code, are there any clear alterations to improve my code?
Hopefully I have explained myself clearly, if not please let me know and I will try to clarify myself. Thank you in advance for the help!
Below is a copy of the code:
function [ascan_output] = TwoD_5_Test1mhz(input_argument)
if(ischar(input_argument))
j = str2double(input_argument);
end
%Above this line is set up for the supercomputer
%%
% create the computational grid
Nx = 2560; % number of grid points in the x (row) direction
Ny = 2560; % number of grid points in the y (column) direction
dx = 7.8125e-6; % grid point spacing in the x direction [m]
dy = 7.8125e-6; % grid point spacing in the y direction [m]
kgrid = kWaveGrid(Nx, dx, Ny, dy);
A=zeros(2560, 2560);
for i = 1:1888
A(i,:)=5000;
end
for i = 1889:2560
A(i,:)=781;
end
for k=j:j+511
A(1892:2560,k)=5000;
end
B=zeros(1280, 2560);
for i = 1:1888
B(i,:)=7300;
end
for i = 1889:2560
B(i,:)=781;
end
for k=j:j+511
B(1892:2560,k)=7300;
end
%%
% define the properties of the propagation medium based on above
medium.sound_speed = A; % [m/s]
medium.density = B; % [kg/m^3]
medium.alpha_coeff = 2.56e-7; % [dB/(MHz^y cm)]
medium.alpha_power = 1.1; %are these right??
t_end = 16e-6; % [s]
kgrid.makeTime(medium.sound_speed, [], t_end);
t_end = 16e-6;
% define the source mask to be a line across the top of the grid
source.p_mask = zeros(Nx, Ny);
source.p_mask(100, 1024:1536) = 1;
%I think the error must be arising from this, as this is the section that has changed?
% define a tone burst and assign it to the x-direction particle velocity
SampleFreq=1/kgrid.dt;
source.p = toneBurst(SampleFreq,1e6, 5, 'Plot', true);
%%
% define a linear sensor along edge of grid
sensor.mask = zeros(Nx, Ny);
sensor.mask(100, 1024:1536) = 1;
sensor.record={'p', 'p_max', 'u'};
% run the simulation with optional inputs for plotting the simulation
% layout in addition to removing the PML from the display
sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, ...
'PMLInside', true, 'PlotLayout', false, 'PlotPML', false,...
'PlotSim', false, 'DataCast', 'single'); %gpuArray-single / single / double
avddatap=mean(sensor_data.p);
timelength=length(avddatap);
time=(1:1:timelength);
time=time*kgrid.dt;
%Below this line is to only save important variables from each simulation
%on supercomputer
clearvars -except avddatap time j
save(sprintf(['TwoD_Test_5_1mhz-',num2str(j)], 'avddatap'));
end