How can I appear the streamline around the spheres?

10 Ansichten (letzte 30 Tage)
Shreen El-Sapa
Shreen El-Sapa am 22 Apr. 2024
Bearbeitet: Voss am 29 Okt. 2024 um 20:12
clc
A=[ -0.2
0.1
-0.4
1.0
-2.1
3.9
-6.9
11.6
-18.7
29.2
-44.0
64.9
-93.5
132.6
-184.2
253.6
-343.1
462.3
-613.4
815.3
-1068.2
1413.3
-1842.6
2462.0
-3232.1
4517.3
-6127.8
10894.9
-17024.9
9869.9];
B=[ 279.2
34.9
-96.3
177.4
-248.7
283.6
-271.0
223.3
-160.5
102.6
-58.7
30.5
-14.4
6.3
-2.5
0.9
-0.3
0.1
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0];
C=[ -4.6
-0.3
0.3
-0.2
0.1
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.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0];
AA=[ -1.4
-0.1
0.4
-1.1
3.0
-7.5
17.7
-40.2
88.1
-188.2
390.9
-797.9
1594.7
-3150.4
6115.4
-11791.2
22393.6
-42435.9
79327.3
-148756.4
275335.2
-515215.3
950981.9
-1800517.1
3351792.8
-6647748.5
12804165.0
-32341874.0
71833152.0
-59216272.0];
BB=[ 15518.5
-406.6
1241.4
-2229.4
3366.8
-4306.2
4732.0
-4611.8
4004.3
-3152.0
2255.2
-1486.4
902.0
-510.4
268.4
-132.9
61.6
-27.2
11.3
-4.5
1.7
-0.6
0.2
-0.1
0.0
0.0
0.0
0.0
0.0
0.0];
CC=[ -15.6
0.3
-0.4
0.3
-0.2
0.1
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.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0];
%%%%%%%%%%%%%%%%%%%%
a = 1 ; %RADIUS
L=.22;
etta=0.01; %0.02;0.05;0.1d0
u2=1; delta=1.5; ettap=.01; alpha=2.d0; %C=a1/a2=0.1
xi=1./sqrt(etta);xi1=1./sqrt(ettap);
alpha1=sqrt((xi.^2+sqrt(xi.^4-4.*xi.^2.*alpha.^2))./sqrt(2));
alpha2=sqrt((xi.^2-sqrt(xi.^4-4.*xi.^2.*alpha.^2))./sqrt(2));
dd=5;
c =-a/L;
b =a/L;
m =a*40; % 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));
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);
%for i1=1:length(x);
% for k1=1:length(x);
% if sqrt(x(i1,k1).^2+y(i1,k1).^2)>1./L;
% r(i1,k1)=0;r2(i1,k1)=0;
% end
% end
%end
warning off
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.^(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,25,'-k','LineWidth',1.1); %,psi2,'--k',psi2,':k'
%[DH1,h1]=contour(x,y,psi1);
%p1=contour(x,y,psi1,[0.3 0.3],'k','LineWidth',1.1); %,'ShowText','on'
%p2=contour(x,y,psi1,[0.4 0.4],'r','LineWidth',1.1);
%p3=contour(x,y,psi1,[0.5 0.5],'g','LineWidth',1.1);
%p4=contour(x,y,psi1,[0.6 0.6],'b','LineWidth',1.1);
%p5=contour(x,y,psi1,[0.7 0.7],'c','LineWidth',1.1);
%p6=contour(x,y,psi1,[0.8 0.8],'m','LineWidth',1.1);
%p7=contour(x,y,psi1,[0.9 0.9],'y','LineWidth',1.1);
p1=contour(x,y,psi1,[-0.001 0.001],'k','LineWidth',1.1); %,'ShowText','on'
p2=contour(x,y,psi1,[-0.005 .005],'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);
%%%%%%%%%%%%%%% $\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')
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 on
%xticklabels([])
%yticklabels([])
%legend('0.01','0.05','0.1','0.4','0.6','0.8','Location','northwest')
%title('$\frac{\beta_1}{a_1\mu}=\frac{a_1\beta_2}{\mu}=1.0,\;R_{H}=1.0,\;\frac{a_2}{a_1}=2.0$','Interpreter','latex','FontSize',12,'FontName','Times New Roman','FontWeight','Normal')
%title('$(a)\;\; R_{H}=1.0,\;\frac{\kappa}{\mu}=4.0$','Interpreter','latex','FontSize',12,'FontName','Times New Roman','FontWeight','Normal')
%%%%%%%%%%%%%%%%%%%%
view([90 90])

Antworten (1)

Gautam
Gautam am 23 Okt. 2024
Hello Shreen,
You can plot the streamlines by using a series of contour functions the way you've been attempting to.
hold on
p11=contour(y,x,psi1,[1 1],'k','LineWidth',1.1);
p21=contour(y,x,psi1,[3 3],'r','LineWidth',1.1);
p31=contour(y,x,psi1,[4 4],'g','LineWidth',1.1);
p41=contour(y,x,psi1,[5 5],'b','LineWidth',1.1);
p51=contour(y,x,psi1,[10 10],'c','LineWidth',1.1);
p61=contour(y,x,psi1,[50 50],'m','LineWidth',1.1);
p71=contour(y,x,psi1,[100 100],'y','LineWidth',1.1);
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(y1,x1,'-k'),'LineWidth',1.1);
fill(y1,x1,'w')
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(y2,x2,'-k'),'LineWidth',1.1);
fill(y2,x2,'w')
axis('equal')
axis off
hold off
  1 Kommentar
shreen elsapa
shreen elsapa am 29 Okt. 2024 um 20:05
Bearbeitet: Voss am 29 Okt. 2024 um 20:12
Can you tell me about this code ?
clc
A =[ -4.56964
0.11715
-0.49082
1.14190
-2.15372
3.63497
-5.67271
8.40030
-11.89370
16.34167
-21.79783
28.55245
-36.60938
46.44017
-57.94571
71.92785
-88.07967
107.83102
-130.44882
158.64197
-190.75777
232.40297
-279.74454
346.19617
-421.83228
549.04956
-694.61816
1157.37988
-1701.72839
928.16284];
B=[ 3.88442
-0.03321
0.02748
-0.00907
0.00189
-0.00029
0.00003
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000];
AA=[ -4.49021
-0.11744
0.56476
-1.16020
2.19876
-3.70012
5.77588
-8.55282
12.10988
-16.63875
22.19414
-29.07158
37.27501
-47.28453
58.99926
-73.23563
89.68111
-109.79158
132.82062
-161.52635
194.22607
-236.62846
284.83081
-352.49063
429.50192
-559.03229
707.24756
-1178.42310
1732.66882
-945.03845];
BB=[ 3.82074
0.03326
-0.03161
0.00921
-0.00193
0.00029
-0.00004
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000
0.00000];
% Parameters
a = 1 ; %RADIUS
L=.1;
alm1=.1; alpha=1; %u2=1;a2=1;beta1=1;beta2=1;delta=3; RH=alpha magnetic; alm1=alpha frequency
alpha1=((alpha^2-alm1^2 * complex(1, 0))^(0.5));
dd = 6;
c =-a/L;
b =a/L;
m =a*200; % NUMBER OF INTERVALS
%[x,y]=meshgrid((c+dd:(b-c)/m:b),(c:(b-c)/m:b)');
[x,y]=meshgrid((c+dd:(b-c)/m:b),(0:(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);AAi=AA(i-1);BBi=BB(i-1);
psi1=psi1+(Ai.*r.^(-i+1)+r.^(1./2).*besselk(i-1./2,r.*alpha1).*Bi ).*gegenbauerC(i,-1./2, cos(t))+(AAi.*r2.^(-i+1)+r2.^(1./2).*besselk(i-1./2,r2.*alpha1).*BBi).*gegenbauerC(i,-1./2,zet);
end
hold on
[DH1,h1]=contour(x,y,psi1,25,'-k','LineWidth',1.1); %,psi2,'--k',psi2,':k'
%[DH1,h1]=contour(x,y,psi1);
%p1=contour(x,y,psi1,[0.3 0.3],'k','LineWidth',1.1); %,'ShowText','on'
%p2=contour(x,y,psi1,[0.4 0.4],'r','LineWidth',1.1);
%p3=contour(x,y,psi1,[0.5 0.5],'g','LineWidth',1.1);
%p4=contour(x,y,psi1,[0.6 0.6],'b','LineWidth',1.1);
%p5=contour(x,y,psi1,[0.7 0.7],'c','LineWidth',1.1);
%p6=contour(x,y,psi1,[0.8 0.8],'m','LineWidth',1.1);
%p7=contour(x,y,psi1,[0.9 0.9],'y','LineWidth',1.1);
%%%%%%%%%%%%%%%
hold on
t3 = linspace(0,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')
hold on
t2 = linspace(0,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')
box on
axis on
xticklabels([])
yticklabels([])
legend('0.01','0.05','0.1','0.4','0.6','0.8','Location','northwest')
Warning: Ignoring extra legend entries.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%55

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Vector Fields 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!

Translated by