Dear,
I have two input vector matrices as variables alf and gam:
alf = [0, 0.03, 0.1, 0.5, 1, 1.5];
gam = [0, 0.03, 0.1, 0.5, 1, 1.5];
Then I have some complicated steps resulting in one value for each possible combination of the two variables alf and gam.
For example alf=0.03 and gam=1 gives one answer. And alf=0.5 and gam=0 gives an other answer etc.
How can I plot this as a 3D plot with a surface with alf and gam as X and Y axes and the outcome Z axis.
if I try:
[gamM, alfM]=meshgrid(gam,alf);
surf(gamM,alfM,log10(T))
it gives me the error saying:
Error using surf (line 71)
Z must be a matrix, not a scalar or vector.
kind regards Nickel

7 件のコメント

Dyuman Joshi
Dyuman Joshi 2021 年 11 月 25 日
編集済み: Dyuman Joshi 2021 年 11 月 25 日
You might be storing T as row/column vector. Change it into a matrix.
Since you mentioned that you get a unique output for each gam and alf, then your T should have 36 elements. Convert it into a 6x6 matrix and try again.
Nickel Blankevoort
Nickel Blankevoort 2021 年 11 月 25 日
Dear Dyuman Joshi,
That is correct I get 36 individual elements for T.
But how do I convert that into a 6x6 matrix?
I run two For loops simultaneously for alf and gam, calculating T.
Kind regards Nickel
Dyuman Joshi
Dyuman Joshi 2021 年 11 月 25 日
for i=1:numel(alf)
for j=1:numel(gam)
T(i,j)= alf(i)*gam(j) %this is just an exmample, use your existing formula
end
end
Nickel Blankevoort
Nickel Blankevoort 2021 年 11 月 25 日
clc;clear
eL=0;
eR=0;
gammaL=-2;
gammaR=-2;
aL=-0.1;
aR=-0.1;
E=linspace(-0,0.0001,1);
kL=acos((eL-E)/(2*gammaL
kR=acos((eR-E)/(2*gammaR));
h=6.582119514*10^-16;
N = 3;
for alf = [0, 0.03, 0.1, 0.5, 1, 1.5];
for gam = [0, 0.03, 0.1, 0.5, 1, 1.5];
H = diag(gam.*ones(1,N-1),-1);
if N>2
H = diag(gam.*ones(1,N-1),-1)+diag(alf.*ones(1,N-2),-2);
else
H = diag(gam.*ones(1,N-1),-1);
end
H = H+H';
a=1;
b=N;
n=1;
g=(((E(n))*eye(length(H))-H))^-1;
gaa=g(a,a);
gbb=g(b,b);
gab=g(a,b);
gba=g(b,a);
vL=(2*gammaL*sin(kL(n)))/h;
vR=(2*gammaR*sin(kR(n)))/h;
gL=(exp(i*kL(n)))/(-gammaL);
gR=(exp(i*kR(n)))/(-gammaR);
sL=aL.^2*gL;
sR=aR.^2*gR;
d=1+sL*gaa+sR*gbb+sR*sL*(gaa*gbb-gab*gba);
t=i*h*sqrt(vL)*aL*gL*(gba/d)*gR*aR*sqrt(vR);
T(:,n)=[abs(t^2)]; %I would like to store the value of this T for every iteration
end
end
Nickel Blankevoort
Nickel Blankevoort 2021 年 11 月 25 日
編集済み: Nickel Blankevoort 2021 年 11 月 25 日
Dear Dyuman,
thank you for the quick responses. I have included my code now for more clearity.
I try to store the values for T for every iteration.
Nickel
Dyuman Joshi
Dyuman Joshi 2021 年 11 月 25 日
I can't see your code? Or did you mean, that you have included what I wrote in your code?
Nickel Blankevoort
Nickel Blankevoort 2021 年 11 月 25 日
Can you see my code now in the attachement?
Nickel

サインインしてコメントする。

 採用された回答

Dyuman Joshi
Dyuman Joshi 2021 年 11 月 25 日

1 投票

I have modified your code -
clc;clear
eL=0;
eR=0;
gammaL=-2;
gammaR=-2;
aL=-0.1;
aR=-0.1;
E=linspace(-0,0.0001,1);
kL=acos((eL-E)/(2*gammaL));
kR=acos((eR-E)/(2*gammaR));
h=6.582119514*10^-16;
N = 3;
alf = [0, 0.03, 0.1, 0.5, 1, 1.5];
gam = [0, 0.03, 0.1, 0.5, 1, 1.5];
for i=1:6
for j=1:6
H = diag(gam(i).*ones(1,N-1),-1);
if N>2
H = diag(gam(i).*ones(1,N-1),-1)+diag(alf(i).*ones(1,N-2),-2);
else
H = diag(gam(i).*ones(1,N-1),-1);
end
H = H+H';
a=1;
b=N;
n=1;
g=(((E(n))*eye(length(H))-H))^-1;
gaa=g(a,a);
gbb=g(b,b);
gab=g(a,b);
gba=g(b,a);
vL=(2*gammaL*sin(kL(n)))/h;
vR=(2*gammaR*sin(kR(n)))/h;
gL=(exp(i*kL(n)))/(-gammaL);
gR=(exp(i*kR(n)))/(-gammaR);
sL=aL.^2*gL;
sR=aR.^2*gR;
d=1+sL*gaa+sR*gbb+sR*sL*(gaa*gbb-gab*gba);
t=i*h*sqrt(vL)*aL*gL*(gba/d)*gR*aR*sqrt(vR);
T(i,j)=[abs(t^2)]; %I would like to store the value of this T for every iteration
end
end
[gamM, alfM]=meshgrid(alf,alf);
surf(gamM,alfM,log10(T))

2 件のコメント

Nickel Blankevoort
Nickel Blankevoort 2021 年 12 月 1 日
thanks a lot for the help, much appriciated!
Nickel Blankevoort
Nickel Blankevoort 2021 年 12 月 1 日
I didn't get the graphs I was expecting to get. Turned out you used i for the for loop too hahah, and in my equations I use i as imaginary... took me a couple of days to figure that out:)
But all happy now!
kind regards and big thanks again

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

ヘルプ センター および File ExchangeMATLAB についてさらに検索

製品

リリース

R2018b

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by