3D plot in polar coordinates

조회 수: 6 (최근 30일)
Hexe
Hexe 2023년 5월 21일
댓글: Star Strider 2023년 5월 23일
Hello! I made a code for solving the integral and it looks realistic in polar coordinates. But, how to present it in 3D as a figure in volume? Perhaps I should add an aditional rotation angle for this. Maybe it is easy for those who know ow to do it. If somebody knows, please, help me. Thank you.
s = 3;
n = 1;
t = 0.1;
r = 1;
a = 0:1:360;
a = a*pi/180;
b = sqrt(2*n*t);
L = sqrt((4*t+r^2)/3);
fun = @(k,u,c,a) ((k.^2).*exp(-1.5*k.^2)).*((u.^2).*(1-u.^2).*exp(-(b.*u).^2).*(cos(s.*k.*u.*cos(a)/L))).*(((cos(c)).^2).*(cos(s.*k.*sqrt(1-u.^2).*sin(a).*(cos(c)+sin(c))/(L*sqrt(2)))));
f3 = arrayfun(@(a)integral3(@(k,u,c)fun(k,u,c,a),0,Inf,-1,1,0,2*pi),a);
B = ((6*sqrt(6)*b^3)/(erf(b)*pi^2))*(1-(3/(2*b^2))*(1-((2*b*exp(-b^2))/(erf(b)*sqrt(pi))))).^(-1);
R = B*f3;
polar(a,R);

채택된 답변

Star Strider
Star Strider 2023년 5월 22일
This took a while to get working correctly, and takes about 500 seconds to run, so I will post the code here slthough not the plot —
s = 3;
n = 1;
tv = 0:0.1:2;%3;
r = 1;
a = 0:2:360;
a = deg2rad(a);
tic
for k = 1:numel(tv)
t = tv(k)
b = sqrt(2*n*t);
L = sqrt((4*t+r^2)/3);
fun = @(k,u,c,a) ((k.^2).*exp(-1.5*k.^2)).*((u.^2).*(1-u.^2).*exp(-(b.*u).^2).*(cos(s.*k.*u.*cos(a)/L))).*(((cos(c)).^2).*(cos(s.*k.*sqrt(1-u.^2).*sin(a).*(cos(c)+sin(c))/(L*sqrt(2)))));
f3 = arrayfun(@(a)integral3(@(k,u,c)fun(k,u,c,a),0,Inf,-1,1,0,2*pi),a);
B = ((6*sqrt(6)*b^3)/(erf(b)*pi^2))*(1-(3/(2*b^2))*(1-((2*b*exp(-b^2))/(erf(b)*sqrt(pi))))).^(-1);
R = B*f3;
ta = t*ones(size(a));
[x,y,z] = pol2cart(a, R, ta);
X(k,:) = x;
Y(k,:) = y;
Z(k,:) = z;
toc % Total Time To This Point
LIT = toc/k % Mean Loop Iteration Time
end
toc
figure
surf(X, Y, Z)
colormap(turbo)
axis('equal')
axis('off')
hold on
xc = [0:0.25*r:r].'*cos(a);
yc = [0:0.25*r:r].'*sin(a);
xr = [0;r]*cos(0:pi/4:2*pi);
yr = [0;r]*sin(0:pi/4:2*pi);
plot3(xc.', yc.', zeros(size(xc)), '-k')
plot3(xr, yr, zeros(size(xr)), '-k')
zt = (0:1:max(tv)).';
zt1 = ones(size(zt));
surf(zt1*xc(end,:), zt1*yc(end,:), zt*ones(size(a)), 'FaceAlpha',0, 'MeshStyle','row')
hold off
text(xr(2,1:end-1)*1.1,yr(2,1:end-1)*1.1, zeros(1,size(xr,2)-1), compose('%3d°',(0:45:315)), 'Horiz','center','Vert','middle')
text(ones(size(zt))*xc(end,end-10), ones(size(zt))*yc(end,end-10), zt, compose('%.1f',zt))
toc
figure
surf(X, Y, Z, 'EdgeColor','none')
colormap(turbo)
axis('equal')
axis('off')
hold on
xc = [0:0.25*r:r].'*cos(a);
yc = [0:0.25*r:r].'*sin(a);
xr = [0;r]*cos(0:pi/4:2*pi);
yr = [0;r]*sin(0:pi/4:2*pi);
plot3(xc.', yc.', zeros(size(xc)), '-k')
plot3(xr, yr, zeros(size(xr)), '-k')
zt = (0:1:max(tv)).';
zt1 = ones(size(zt));
surf(zt1*xc(end,:), zt1*yc(end,:), zt*ones(size(a)), 'FaceAlpha',0, 'MeshStyle','row')
hold off
text(xr(2,1:end-1)*1.1,yr(2,1:end-1)*1.1, zeros(1,size(xr,2)-1), compose('%3d°',(0:45:315)), 'Horiz','center','Vert','middle')
text(ones(size(zt))*xc(end,end-10), ones(size(zt))*yc(end,end-10), zt, compose('%.1f',zt))
toc
TTmin = toc/60
% END
I only run it here from 0 to 2, because there is not much detail beyond that. The time dimension is the ‘Z’ dimension. Most of the detail is between 0 and 0.5.
.
  댓글 수: 6
Hexe
Hexe 2023년 5월 23일
편집: Hexe 2023년 5월 23일
Dear Star Strider.
Your previous answer is good and acceptable. I just wanted to clarify about another way of presentation the result, but probably it is not possible.
Thank you a lot.
Star Strider
Star Strider 2023년 5월 23일
As always, my pleasure!

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Graphics에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by