I want to plot the half only in the attached figure. How can I do that?
5 ビュー (過去 30 日間)
古いコメントを表示
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%subplot(1,2,1)
A =[ -16.
0.
0.
0.
0.
1.
-2.
7.
-20.
54.
-145.
384.
-989.
2525.
-6332.
15811.
-38887.
95760.
-232600.
569166.
-1374036.
3371349.
-8148474.
20348640.
-49802872.
131423232.
-334120288.
1149872256.
-3565926912.
4013070592.];
B=[ -350738.
26.
-101.
194.
-291.
368.
-403.
395.
-347.
280.
-207.
142.
-91.
54.
-31.
16.
-8.
4.
-2.
1.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.];
C=[ 35.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.];
AA=[ -4.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
-1.
1.
-1.
1.
-1.
1.
-2.
2.
-3.
3.
-4.
5.
-9.
14.
-8.];
BB=[ -95.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.];
CC=[ 5.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.
0.];
a = 1 ; %RADIUS
L=.1;
akm=4;gamma=0.3;arh=1; %beta1=beta2=1,a1=1,a2=2,arh=1,delta=0.1,u2=1
alphaa=sqrt(((2+akm).*akm./(gamma.*(2+akm))).^2+arh.^2);
betaa=(2.*akm.*arh.^2./gamma).^(0.25);
alpha1=sqrt((alphaa.^2+sqrt(alphaa.^4-4.*betaa.^4))./2);
alpha2=sqrt((alphaa.^2-sqrt(alphaa.^4-4.*betaa.^4))./2);
dd=5;
c =-a/L;
b =a/L;
m =a*60; % NUMBER OF INTERVALS
[x,y]=meshgrid((c+dd:(b-c)/m:b),(c:(b-c)/m:b)');
[I, J]=find(sqrt(x.^2+y.^2)<(a-0.1));
if ~isempty(I)
x(I,J) = 0; y(I,J) = 0;
end
r=sqrt(x.^2+y.^2);
t=atan2(y,x);
r2=sqrt(r.^2+dd.^2-2.*r.*dd.*cos(t));
zet=(r.^2-r2.^2-dd.^2)./(2.*r2.*dd);
warning on
psi1=0;
for i=2:7
Ai=A(i-1);Bi=B(i-1);Ci=C(i-1);AAi=AA(i-1);BBi=BB(i-1);CCi=CC(i-1);
%psi1=-psi1-(Ai.*r.^(-i-1)+r.^(-3./2).*besselk(i-1./2,r.*alpha1).*Bi+r.^(-3./2).*besselk(i-1./2,r.*alpha2).*Ci).*legendreP(i-1,cos(t))-(AAi.*r2.^(-i-1)+r2.^(-3./2).*besselk(i-1./2,r2.*alpha1).*BBi+r2.^(1./2).*besselk(i-1./2,r2.*alpha2).*CCi).*legendreP(i-1,zet);
psi1=psi1+(Ai.*r.^(-i+1)+r.^(1./2).*besselk(i-1./2,r.*alpha1).*Bi+r.^(1./2).*besselk(i-1./2,r.*alpha2).*Ci).*gegenbauerC(i,-1./2, cos(t))+(AAi.*r2.^(-i+1)+r2.^(1./2).*besselk(i-1./2,r2.*alpha1).*BBi+r2.^(1./2).*besselk(i-1./2,r2.*alpha2).*CCi).*gegenbauerC(i,-1./2,zet);
end
hold on
%[DH1,h1]=contour(x,y,psi1,15,'-k'); %,psi2,'--k',psi2,':k'
p11=contour(x,y,psi1,[1 1],'k','LineWidth',1.1);
p21=contour(x,y,psi1,[3 3],'r','LineWidth',1.1);
p31=contour(x,y,psi1,[4 4],'g','LineWidth',1.1);
p41=contour(x,y,psi1,[5 5],'b','LineWidth',1.1);
p51=contour(x,y,psi1,[10 10],'c','LineWidth',1.1);
p61=contour(x,y,psi1,[50 50],'m','LineWidth',1.1);
p71=contour(x,y,psi1,[100 100],'y','LineWidth',1.1);
%axis square;
%beta1=beta2=1,a1=1,a2=2,arh=1,delta=0.1,u2=1
%akm=4;
%title('$\frac{\beta_1}{a_1\mu}=\frac{a_1\beta_2}{\mu}=1.0,\;R_{H}=1.0,\;\frac{a_2}{a_1}=0.5$','Interpreter','latex','FontSize',12,'FontName','Times New Roman','FontWeight','Normal')
%title('$(b)$','Interpreter','latex','FontSize',12,'FontName','Times New Roman','FontWeight','Normal')
%%%%%%%%%%%%%%% $\frac{\textstyle a_1+a_2}{\textstyle h}=6.0,\;
hold on
t3 = linspace(0,2*pi,1000);
h2=0;
k2=0;
rr2=2;
x2 = rr2*cos(t3)+h2;
y2 = rr2*sin(t3)+k2;
set(plot(x2,y2,'-k'),'LineWidth',1.1);
fill(x2,y2,'w')
%axis square;
hold on
t2 = linspace(0,2*pi,1000);
h=dd;
k=0;
rr=1;
x1 = rr*cos(t2)+h;
y1 = rr*sin(t2)+k;
set(plot(x1,y1,'-k'),'LineWidth',1.1);
fill(x1,y1,'w')
%axis square;
axis('equal')
axis off
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 件のコメント
回答 (1 件)
Gautam
2024 年 10 月 23 日 5:42
Hello @Shreen El-Sapa
You can make the following changes to the contour plot to plot only the top half
p1=contour(x,y,psi1,[0.01 0.01],'k','LineWidth',1.1); %,'ShowText','on'
p2=contour(x,y,psi1,[0.05 .05],'r','LineWidth',1.1);
p3=contour(x,y,psi1,[0.1 0.1],'g','LineWidth',1.1);
p4=contour(x,y,psi1,[0.4 0.4],'b','LineWidth',1.1);
p5=contour(x,y,psi1,[0.6 0.6],'c','LineWidth',1.1);
p6=contour(x,y,psi1,[0.8 0.8],'m','LineWidth',1.1);
And change the variables "t2" and "t3" to:
t2 = linspace(0,pi,1000);
t3 = linspace(0,2*pi,1000);
0 件のコメント
参考
カテゴリ
Help Center および File Exchange で Axis Labels についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!