Help with a contour graph
2 views (last 30 days)
Show older comments
I've been trying (and failing to plot) this graph (Fig 1) for a week, but the closest I've gotten is Fig 2. It seems weirdly stretched out and I can't figure out what's wrong.
This is my code:
clear all
close all
clc
%Constant
rho = 4420; %kg/m^3
Cp = 550; %J/kg?K
T0 = 303.15; %K
A = 0.5; %[Absorbtivity]
k = 7.2; %W/m/K
alpha = 2.96*10^-6; %m^2/s
D = alpha;
P = 100; %W
v = 1; %m/s
u = v;
Tm = 1933; %K
d_laser = 0.01; %mm
r_laser = d_laser/2; %mm
a = r_laser;
p = D/(u*a);
%Define
x = linspace(-0.05,0.125,100);
y = linspace(-0.025,0.025,100);
z = linspace(0,0.05,100);
%Normalize
x_nor = x/a;
y_nor = y/a;
z_nor = z/sqrt((D*a/u))
[x_mesh,y_mesh,z_mesh] = ndgrid(x_nor,y_nor,z_nor);
% fun = @(t) exp((-z_mesh.^2/(4*t)-((y_mesh.^2.+(x_mesh-t).^2)./(4*p.*t+1)))./((4*p*t)+1)*sqrt(t));
fun = @(t) exp((-z_mesh.^2/(4*t))-((y_mesh.^2+(x_mesh-t).^2)/((4*p*t)+1)))/(((4*p*t)+1)*sqrt(t));
g = integral(fun,0,Inf,'ArrayValued',true);
squeeze(g(:,:,1))
figure(1)
iz = 1;
contourf(x_mesh(:,:,iz), y_mesh(:,:,iz),squeeze(g(:,:,iz)))
axis([-5 25 -5 5])
% xlim([5 25])
% ylim([-10 0])
colorbar
% figure(2)
% iy =1;
% contourf(x_mesh,z_mesh,transpose(squeeze(g(:,iy,:))))
% axis([-5 25 -10 0])
% colorbar
4 Comments
Sargondjani
on 29 Jun 2023
You can set the levels at which a contour should appear. Did you try setting that?
(Sorry, it was too much code for me to fnd it quickly)
cdawg
on 29 Jun 2023
The paper says 0.1 mm and it looks like you have D = 0.01 mm. Not sure if this helps
Answers (1)
akshatsood
on 22 Aug 2023
Hi Sarinya,
I thoroughly investigated the code script attached and have concluded that the weirdly stretch visible in the graph is due to improper scaling of the dimensions. For instance, consider the following line from the code snippet.
p = D/(u*a)
D has units m^2/s, u has units m/s but a has units mm. So, for the expression to be evaluated in a correct manner, uniformity in dimensions has to be maintained. To correct this, a simple modification shown below can be performed.
p = 1000*D/(u*a)
On similar lines, a modification while computing z_nor needs to be done
z_nor = z/sqrt((1000*D*a/u))
I also went through the research paper attached and noticed that diameter of the laser is stated to be 0.1 mm as highlighted in the comments as well. Performing the following changes would be sufficient to replicate the required results.
d_laser = 0.1; %mm
%Define
x = linspace(-0.5,1.25,100);
y = linspace(-0.25,0.25,100);
z = linspace(0,0.5,100);
I hope this helps.
0 Comments
See Also
Categories
Find more on 2-D and 3-D Plots in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!