hi dear professor bradley treeby
I considered the heat source as a 3-d Gaussian function and I want to know the temperature in the xy-plane in front of HIFU.but I think the temperature after heating and cooling is unreasonable.may you tell me what is the problem?
thanks in advance
here is my code
%% set the size of the perfectly matched layer (PML)
PML_X_SIZE = 10; % [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 = 120- 2*PML_X_SIZE; % [grid points]
Ny = 120- 2*PML_Y_SIZE; % [grid points]
Nz = 120- 2*PML_Z_SIZE; % [grid points]
%% set desired grid size in the x-direction not including the PML
x = .1; % [m]
%% calculate the spacing between the grid points
dx = x/Nx; % [m]
dy = dx; % [m]
dz = dx; % [m]
%% create the k-space grid
kgrid =kWaveGrid(Nx, dx, Ny, dy, Nz, dz);
x1=linspace(0,.1,100);
y1=x1;
z1=x1;
[X,Y,Z] = meshgrid(x1,y1,z1);
Q=((2.34e4)*exp(-power(((X-.05)/.004284),2))).*((2.34e4)*exp(-power(((Y-.05)/.004284),2))).*((2.008e4)*exp(-power(((Z-.01153)/.02147),2)));
clear medium source sensor;
% set the background temperature and heating term
source.Q = Q;
source.T0 = 37;
% define medium properties related to diffusion
medium.density = 1020; % [kg/m^3]
medium.thermal_conductivity = 0.5; % [W/(m.K)]
medium.specific_heat = 3600; % [J/(kg.K)]
% create kWaveDiffusion object
kdiff = kWaveDiffusion(kgrid, medium, source, []);
% set source on time and off time
on_time = 1; % [s]
off_time = 20; % [s]
% set time step size
dt = 0.1;
% take time steps
kdiff.takeTimeStep(round(on_time / dt), dt);
% store the current temperature field
T1 = kdiff.T;
% turn off heat source and take time steps
kdiff.Q = 0;
kdiff.takeTimeStep(round(off_time / dt), dt);
% store the current temperature field
T2 = kdiff.T;
% =========================================================================
% VISUALISATION
% =========================================================================
% plot the thermal dose and lesion map
figure;
% plot the volume rate of heat deposition
subplot(3, 1, 1);
imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, Q(:,:,50) * 1e-7);
h = colorbar;
xlabel(h, '[kW/cm^3]');
ylabel('x-position [mm]');
xlabel('y-position [mm]');
axis image;
title('Volume Rate Of Heat Deposition');
% plot the temperature after heating
subplot(3, 1, 2);
imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, T1(:,:,50));
h = colorbar;
xlabel(h, '[degC]');
ylabel('x-position [mm]');
xlabel('y-position [mm]');
axis image;
title('Temperature After Heating');
% plot the temperature after cooling
subplot(3, 1, 3);
imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, T2(:,:,50));
h = colorbar;
xlabel(h, '[degC]');
ylabel('x-position [mm]');
xlabel('y-position [mm]');
axis image;
title('Temperature After Cooling');