How can I correct this code?

Question

0 votes

% Data Initialization
A = [ -0.1880
    0.0100
    -0.0079
    -0.0273
    0.1587
    -0.4069
    0.8943
    -1.7404
    3.0826
    -5.1443
    8.1501
    -12.4334
    18.3324
    -26.3551
    36.9611
    -50.9219
    68.8245
    -91.8411
    120.6744
    -157.1298
    201.9740
    -258.0539
    326.0616
    -410.6317
    512.0566
    -638.1943
    788.2017
    -976.1147
    1198.2528
    -1481.4827
    1815.1080
    -2255.8411
    2774.5508
    -3510.8018
    4379.5371
    -5837.7578
    7571.4097
    -12936.2158
    19555.4707
    -10974.3643];
B = [  11665.0576
    -4.7293
    13.9157
    -26.5743
    37.1578
    -42.1387
    44.1123
    -45.5787
    44.7030
    -39.4265
    30.5014
    -20.7728
    12.5115
    -6.7371
    3.2495
    -1.4121
    0.5494
    -0.1905
    0.0572
    -0.0141
    0.0023
    0.0001
    -0.0003
    0.0002
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
C = [   -22.4193
    -0.0217
    0.0677
    -0.1272
    0.1340
    -0.0548
    -0.0433
    0.1017
    -0.1076
    0.0619
    -0.0047
    -0.0195
    0.0135
    -0.0030
    -0.0010
    0.0009
    -0.0002
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
D = [    -0.0871
    -0.0006
    0.0002
    0.0004
    -0.0012
    0.0016
    -0.0018
    0.0017
    -0.0015
    0.0013
    -0.0010
    0.0008
    -0.0006
    0.0004
    -0.0003
    0.0002
    -0.0001
    0.0001
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
E = [    -0.1880
    0.0100
    -0.0079
    -0.0273
    0.1587
    -0.4069
    0.8943
    -1.7404
    3.0826
    -5.1443
    8.1501
    -12.4334
    18.3324
    -26.3551
    36.9611
    -50.9219
    68.8245
    -91.8411
    120.6744
    -157.1298
    201.9740
    -258.0539
    326.0616
    -410.6317
    512.0566
    -638.1943
    788.2017
    -976.1147
    1198.2528
    -1481.4827
    1815.1080
    -2255.8411
    2774.5508
    -3510.8018
    4379.5371
    -5837.7578
    7571.4097
    -12936.2158
    19555.4707
    -10974.3643];
F = [  11665.0576
    -4.7293
    13.9157
    -26.5743
    37.1578
    -42.1387
    44.1123
    -45.5787
    44.7030
    -39.4265
    30.5014
    -20.7728
    12.5115
    -6.7371
    3.2495
    -1.4121
    0.5494
    -0.1905
    0.0572
    -0.0141
    0.0023
    0.0001
    -0.0003
    0.0002
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
% Parameters
a = 1;                         % Radius
L = 0.1;                       % Length
dd = 6;                        % Separation distance
alpha = 10;
etta = 0.01;
ettap = 0.1;                   % Additional parameters
zalm = alpha * complex(1, 0) / sqrt(2.d0);
xi = 1.0 / etta;
xi1 = 1.0 / ettap;
zk1 = sqrt((xi + sqrt(xi^2 - 4.d0 * xi * zalm^2)) / 2.0);
zk2 = sqrt((xi - sqrt(xi^2 - 4.d0 * xi * zalm^2)) / 2.0);
c = -a / L;
b = a / L;
m = a * 50;                    % Number of intervals
% Create Mesh Grid
[x, y] = meshgrid(c + dd : (b - c) / m : b, c : (b - c) / m : b);
% Apply Boundary Conditions
[I, J] = find(sqrt(x.^2 + y.^2) < (a - 0.01));
if ~isempty(I)
    x(I, J) = NaN; % Use NaN to avoid affecting the visualization
    y(I, J) = NaN;
end
% Compute Values
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);
% Calculate Psi1
psi1 = zeros(size(x));
for i = 2:5
    Ai = A(i-1);
    Bi = B(i-1);
    Ci = C(i-1);
    Di = D(i-1);
    Ei = E(i-1);
    Fi = F(i-1);
    psi1 = psi1 + (Ai .* r.^(-i + 1) + r.^(1 / 2) .* besselk(i - 1 / 2, r .* zk1) .* Bi + ...
        r.^(1 / 2) .* besselk(i - 1 / 2, r .* zk2) .* Ci) .* ...
        gegenbauerC(i, -1 / 2, cos(t)) + ...
        (Di .* r2.^(-i + 1) + r2.^(1 / 2) .* besselk(i - 1 / 2, r2 .* zk1) .* Ei + ...
        r2.^(1 / 2) .* besselk(i - 1 / 2, r2 .* zk2) .* Fi) .* ...
        gegenbauerC(i, -1 / 2, zet);
end
% Plot Contours with Different Colors
figure;
contourf(y, x, psi1, 20); % Create filled contours with automatic color mapping
Error using contourf (line 55)
Input arguments must be real. Use the function REAL to get the real part of the inputs.
colorbar; % Add colorbar to check values
axis equal;
title('$\delta=1.5,\frac{\beta\mu}{a_1}=1.0,\alpha=1.0,\ell=0.1$', ...
    'Interpreter', 'latex', 'FontSize', 12, 'FontName', 'Times New Roman', 'FontWeight', 'Normal');
% Plot Spheres
hold on;
% Sphere 1 (Vertical)
t3 = linspace(0, 2 * pi, 1000);
rr2 = 2.5;
x2 = rr2 * cos(t3); % Keep x and y aligned for vertical spheres
y2 = rr2 * sin(t3);
plot(y2, x2, '-k', 'LineWidth', 1.1);
fill(y2, x2, [0.5 0.5 0.5]); % Fill with gray color
% Sphere 2 (Vertical)
t2 = linspace(0, 2 * pi, 1000);
h = dd;
rr = 2.5;
x1 = rr * cos(t2) + h; % Keep x and y aligned for vertical spheres
y1 = rr * sin(t2);
plot(y1, x1, '-k', 'LineWidth', 1.1);
fill(y1, x1, [0.5 0.5 0.5]); % Fill with gray color
% Axis Formatting
axis off;

7 Comments
Show 5 older comments Hide 5 older comments

Shreen El-Sapa on 20 Mar 2025

Edited: Cris LaPierre on 20 Mar 2025

Open in MATLAB Online

I changed it to be real but no contours appear.

A  =[       -0.3959
    -0.0099
    0.1405
    -1.6836
    5.0520
    -13.4187
    29.7987
    -58.8179
    107.2008
    -184.0086
    299.9961
    -469.8972
    709.5017
    -1041.9386
    1489.3081
    -2087.1648
    2864.5986
    -3876.0586
    5157.6006
    -6793.5161
    8825.0547
    -11385.5527
    14515.8311
    -18433.4863
    23164.8926
    -29079.9844
    36157.7695
    -45061.1953
    55644.1328
    -69180.5000
    85205.2031
    -106419.5000
    131503.4375
    -167137.9531
    209373.5938
    -280203.1563
    364797.9688
    -625535.8750
    948878.1875
    -534257.8750];
B=[       0.9706
    -0.0808
    0.1285
    0.0734
    -0.3791
    0.5158
    -0.4338
    0.2585
    -0.1138
    0.0378
    -0.0093
    0.0014
    0.0000
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
C=[        1.5865
    -0.0136
    0.0671
    -0.0500
    0.1636
    -0.1284
    -0.0738
    0.1030
    -0.0130
    -0.0151
    0.0051
    0.0003
    -0.0004
    0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
AA=[     -0.0495
    0.0018
    -0.0054
    0.0274
    -0.0411
    0.0547
    -0.0610
    0.0605
    -0.0553
    0.0474
    -0.0386
    0.0302
    -0.0228
    0.0167
    -0.0119
    0.0084
    -0.0057
    0.0039
    -0.0026
    0.0017
    -0.0011
    0.0007
    -0.0005
    0.0003
    -0.0002
    0.0001
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
BB=[     -0.3959
    -0.0099
    0.1405
    -1.6836
    5.0520
    -13.4187
    29.7987
    -58.8179
    107.2008
    -184.0086
    299.9961
    -469.8972
    709.5017
    -1041.9386
    1489.3081
    -2087.1648
    2864.5986
    -3876.0586
    5157.6006
    -6793.5161
    8825.0547
    -11385.5527
    14515.8311
    -18433.4863
    23164.8926
    -29079.9844
    36157.7695
    -45061.1953
    55644.1328
    -69180.5000
    85205.2031
    -106419.5000
    131503.4375
    -167137.9531
    209373.5938
    -280203.1563
    364797.9688
    -625535.8750
    948878.1875
    -534257.8750];
CC=[      0.9706
    -0.0808
    0.1285
    0.0734
    -0.3791
    0.5158
    -0.4338
    0.2585
    -0.1138
    0.0378
    -0.0093
    0.0014
    0.0000
    -0.0001
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000
    0.0000];
a = 1 ;   %RADIUS
L=.11;
alpha=10; etta=0.3; ettap=0.1;     %alpha=lambda; alpha1=alpha
lambda=real(complex(0,alpha))/sqrt(2); 
xi = 1.0 / etta;
xi1 = 1.0 / ettap;
alpha1 = sqrt((xi + sqrt(xi^2 - 4.d0 * xi * lambda^2)) / 2.0);
alpha2 = sqrt((xi - sqrt(xi^2 - 4.d0 * xi * lambda^2)) / 2.0);
c =-a/L;
b =a/L;
m =a*50; % NUMBER OF INTERVALS
[x,y]=meshgrid([c+dd:(b-c)/m:b],[c:(b-c)/m:b]');
Unrecognized function or variable 'dd'.
[I J]=find(sqrt(x.^2+y.^2)<(a-.01));
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 off
psi1=0;
for i=2:10;
    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,real(psi1),20,'-k');
%axis square;
axis equal
axis equal
%%%%%%%%%%%%%%%%  $\frac{\textstyle a_1+a_2}{\textstyle h}=6.0,\;
hold on
t3 = linspace(0,2*pi,1000);
h2=0;
k2=0;
rr2=1;
x2 = rr2*cos(t3)+h2;                   
y2 = rr2*sin(t3)+k2;
set(plot(x2,y2,'-k'),'LineWidth',1.1);
fill(x2,y2,'w')
%axis square;
axis equal
hold on
t2 = linspace(0,2*pi,1000);
h=dd;
k=0;
rr=2;
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
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Shreen El-Sapa on 20 Mar 2025

Edited: Walter Roberson on 20 Mar 2025

Open in MATLAB Online

clc

A =[ -0.3959

-0.0099

0.1405

-1.6836

5.0520

-13.4187

29.7987

-58.8179

107.2008

-184.0086

299.9961

-469.8972

709.5017

-1041.9386

1489.3081

-2087.1648

2864.5986

-3876.0586

5157.6006

-6793.5161

8825.0547

-11385.5527

14515.8311

-18433.4863

23164.8926

-29079.9844

36157.7695

-45061.1953

55644.1328

-69180.5000

85205.2031

-106419.5000

131503.4375

-167137.9531

209373.5938

-280203.1563

364797.9688

-625535.8750

948878.1875

-534257.8750];

B=[ 0.9706

-0.0808

0.1285

0.0734

-0.3791

0.5158

-0.4338

0.2585

-0.1138

0.0378

-0.0093

0.0014

0.0000

-0.0001

0.0000

0.0000];

C=[ 1.5865

-0.0136

0.0671

-0.0500

0.1636

-0.1284

-0.0738

0.1030

-0.0130

-0.0151

0.0051

0.0003

-0.0004

0.0001

0.0000

0.0000];

AA=[ -0.0495

0.0018

-0.0054

0.0274

-0.0411

0.0547

-0.0610

0.0605

-0.0553

0.0474

-0.0386

0.0302

-0.0228

0.0167

-0.0119

0.0084

-0.0057

0.0039

-0.0026

0.0017

-0.0011

0.0007

-0.0005

0.0003

-0.0002

0.0001

-0.0001

0.0000

0.0000];

BB=[ -0.2125

-0.0273

0.1340

-0.1003

0.0364

-0.0009

-0.0062

0.0034

-0.0011

0.0002

0.0000

0.0000];

CC=[ -0.1629

0.0227

0.0273

-0.0859

0.0640

-0.0009

-0.0110

0.0027

0.0002

-0.0001

0.0000

0.0000];

a = 1 ; %RADIUS

L=.11;

etta=0.3; ettap=0.1; alpha=10;delta=1.5;

lambda= complex(0,alpha)/sqrt(2);

xi=sqrt(1./etta); xi1=sqrt(1./ettap);

alpha1=sqrt((xi+sqrt(xi^2-4.*xi.*lambda.^2))./2);

alpha2=sqrt((xi-sqrt(xi^2-4.*xi.*lambda.^2))./2);

dd=6; %h

c =-a/L;

b =a/L;

m =a*50; % NUMBER OF INTERVALS

%[x,y]=meshgrid((c+dd:(b-c)/m:b),(c:(b-c)/m:b)');

[x,y]=meshgrid((c:(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) = NaN; % Use NaN to avoid affecting the visualization

y(I, J) =NaN;

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:15

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,real(psi1),'-k','LineWidth',1.1); %,psi2,'--k',psi2,':k'

%[DH1,h1]=contour(x,y,psi1);

%%%%%%%%%%%%%%%

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,[0.7 0.7 0.7])

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,[0.7 0.7 0.7])

axis off

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Shreen El-Sapa on 20 Mar 2025

Edited: Walter Roberson on 21 Mar 2025

Open in MATLAB Online

clc

A =[ 0.4383

-0.4080

1.4130

-3.7497

8.7052

-18.3631

35.6095

-64.5347

110.4802

-180.8033

284.0856

-432.3704

638.2991

-920.6697

1296.9564

-1795.9568

2440.4854

-3274.5806

4326.2554

-5663.6147

7318.2793

-9397.8857

11932.9268

-15098.8994

18913.5957

-23675.0488

29361.4824

-36506.2852

44985.0234

-55821.1367

68631.0156

-85581.2891

105597.8906

-134030.5938

167689.9219

-224158.2500

291518.9688

-499382.0625

756808.3125

-425743.2188];

B=[ 784.3210

-41.3373

140.2780

-272.2328

368.2668

-379.0945

310.1994

-209.5842

118.3485

-54.6368

18.1556

-1.5761

-3.6361

3.8433

-2.6156

1.4482

-0.6952

0.2997

-0.1176

0.0427

-0.0144

0.0045

-0.0013

0.0004

-0.0001

0.0000

0.0000];

C=[ 3.9761

-0.1226

0.3397

-0.5729

0.5937

-0.3249

0.0574

0.0449

-0.0379

0.0126

-0.0014

-0.0005

0.0002

0.0000

0.0000];

AA=[ -0.1536

0.0249

-0.0437

0.0587

-0.0688

0.0729

-0.0710

0.0645

-0.0553

0.0453

-0.0356

0.0271

-0.0200

0.0144

-0.0101

0.0070

-0.0048

0.0032

-0.0021

0.0014

-0.0009

0.0006

-0.0004

0.0002

-0.0001

0.0001

-0.0001

0.0000

0.0000];

BB=[ 6.3555

-0.2784

-0.1038

1.4075

-2.6351

2.8560

-2.3294

1.5033

-0.8018

0.3637

-0.1427

0.0493

-0.0152

0.0042

-0.0011

0.0002

-0.0001

0.0000

0.0000];

CC=[ -1.2158

0.1330

-0.2657

0.1976

-0.0445

-0.0175

0.0124

-0.0028

0.0001

0.0000

0.0000E+00

0.0000E+00];

a = 1 ; %RADIUS

L=.11;

etta=0.01; ettap=0.01; alpha=4;delta=1.5;

lambda= complex(0,alpha)/sqrt(2);

xi=sqrt(1./etta); xi1=sqrt(1./ettap);

alpha1=sqrt((xi+sqrt(xi^2-4.*xi.*lambda.^2))./2);

alpha2=sqrt((xi-sqrt(xi^2-4.*xi.*lambda.^2))./2);

dd=6; %h

c =-a/L;

b =a/L;

m =a*50; % 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.01));

if ~isempty(I)

x(I,J) = NaN; % Use NaN to avoid affecting the visualization

y(I, J) =NaN;

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:25

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,real(psi1),'--k','LineWidth',1.1,'ShowText','on'); %,psi2,'--k',psi2,':k'

[DH1,h1]=contour(y,x,real(psi1),'-k');

%p1=contour(y,x,real(psi1),[0.001 0.001],'-c','LineWidth',1.3,'ShowText','on'); %,'ShowText','on'

%p2=contour(y,x,real(psi1),[0.005 0.005],'-r','LineWidth',1.3);

%p3=contour(y,x,real(psi1),[0.05 0.05],'-b','LineWidth',1.3);

%p4=contour(y,x,real(psi1),[.1 .1],'-k','LineWidth',1.3);

%p5=contour(y,x,real(psi1),[.2 .2],'--c','LineWidth',1.3);

%p6=contour(y,x,real(psi1),[.4 .4],'--r','LineWidth',1.3);

%p7=contour(y,x,real(psi1),[.6 .6],'--b','LineWidth',1.3);

%p8=contour(y,x,real(psi1),[.8 .8],'--k','LineWidth',1.3);

%p9=contour(y,x,real(psi1),[1 1],'-.c','LineWidth',1.3);

%p10=contour(y,x,real(psi1),[2 2],'-.r','LineWidth',1.3);

%p11=contour(y,x,real(psi1),[3 3],'-.b','LineWidth',1.3);

%p12=contour(y,x,real(psi1),[4 4],'-.k','LineWidth',1.3);

%axis equal; % Maintain aspect ratio

%axis tight; % Adjust axis limits to fit the data

%%%%%%%%%%%%%%%

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,[0.7 0.7 0.7])

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,[0.7 0.7 0.7])

axis off

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Torsten on 21 Mar 2025

Open in MATLAB Online

Try

[DH1,h1]=contour(y,x,real(psi1));
colorbar

and you will see that the range of your values for psi1 goes up until 1e154.

So a variation of the contour colors is restricted to a very small region in the graphics.

Sign in to comment.

Sign in to answer this question.

Follow Question

How can I correct this code?

7 Comments
Show 5 older comments Hide 5 older comments

Answers (0)

Categories

Tags

Community Treasure Hunt

How can I correct this code?

7 Comments Show 5 older comments Hide 5 older comments

Answers (0)

Categories

Tags

See Also

Community Treasure Hunt

7 Comments
Show 5 older comments Hide 5 older comments