I want to plot the half only in the attached figure. How can I do that?
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%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 Kommentare
Antworten (1)
Gautam
am 23 Okt. 2024
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 Kommentare
Siehe auch
Kategorien
Mehr zu Line Plots finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!