Please find the mistake in this code.
Show older comments
Reference image: plot should be like this but don't know where I am making mistake.syms x lambda
warning off
alpha = -0.1;
sigma = 0.1;
eps = -0.1;
e = 0.2;
M = 10;
a = 2;
figure
psi_list = [10, 20, 30, 40. 50];
for i = 1:numel(psi_list)
psi = psi_list(i);
hbar = @(x) a - a.*x + x;
A1 = eps + alpha^3 + (3 * sigma^2 * alpha);
B1 =@(x,lambda) (-3 * lambda * M) * ((hbar(x).^2) + (2 .* hbar(x) .* alpha) + (sigma^2) + (alpha^2));
a1 =@(x) tanh(M .* hbar(x));
b1 =@(x) 1 - ((tanh(M .* hbar(x))).^2);
c1 = (M * alpha) - ((M^3 * A1)/3);
d1 = @(lambda) 2 * (M^2) * (1 + lambda);
C1=@(x) a1(x) + (b1(x) .* c1);
D1 =@(x,lambda) d1(lambda) .* ((hbar(x).^3) + (3 .* (hbar(x).^2) .* alpha) + (3 .* hbar(x) .* (alpha)^2) + (3 .* hbar(x) .* (sigma)^2) + eps + (3 * alpha * (sigma^2)) + (alpha^3));
f1 =@(x,lambda) B1(x,lambda) + (3 * lambda .* C1(x) .* hbar(x)) + (3 * lambda .* C1(x) .* alpha) + (D1(x,lambda) .* C1(x));
f2 =@(x,lambda) 12 * (M^2) * (1 + lambda) .* C1(x);
f3 = psi * (e^3);
f4 = @(lambda) (1 + lambda) *180 * ((1 - e)^2);
f5 =@(lambda) 1/(2 + lambda);
F = @(x,lambda) ((f5(lambda) .* f1(x,lambda))./f2(x,lambda)) + (f3/f4(lambda));
q1 =@(x,lambda) hbar(x) ./ (2 .* F(x,lambda));
Q1 = int(q1,x,0,1); % integration with respect to "x" limit 0 to 1
q2 = @(x,lambda) 1./(F(x,lambda));
Q2 = int(q2,x,0,1); % integration with respect to "x" limit 0 to 1
Q = @(lambda) Q1(lambda)/Q2(lambda); % Q is a function of "lambda"
p1 = @(x,lambda) (1./F(x,lambda)) .* ((0.5 .* hbar(x)) - Q(lambda));
P = @(x,lambda) int(p1,x); % integration with respect to "x" limit 0 to x
% W is obtained by integrating "P" with respect to x limit to x is 0 to
% 1
W =@(lambda) int(P,x,0,1);
fplot(@(lambda) W(lambda))
end
legend(num2str(psi_list'))
ylim([0 1])
xlim([0 4])
set(gca, 'ytick', 0:0.1:1);
set(gca, 'xtick', 0:0.2:4);
xlabel('lambda')
ylabel('W(lambda)')
6 Comments
dpb
on 8 Jul 2022
No have symbolic TB here, sorry.
I'd suggest working through logic w/ debugger to see where it fails to do as expected first...
Steven Lord
on 8 Jul 2022
warning off
I strongly recommend you not do this. MATLAB doesn't issue warnings "just because", it issues warnings because it sees something in the code that may cause a problem later or may indicate that the code isn't doing what you think it's doing.
So when you eliminate that line of code, what warnings get issued? Show the full and exact text of those messages (all the text displayed in orange in the Command Window.) Those warnings may very well explain why you're not seeing anything on your axes.
AVINASH SAHU
on 8 Jul 2022
syms x lambda
% warning off
alpha = -0.1;
sigma = 0.1;
eps = -0.1;
e = 0.2;
M = 10;
a = 2;
figure
hold on
psi_list = [10, 20, 30, 40. 50];
for i = 1:numel(psi_list)
psi = psi_list(i);
hbar = @(x) a - a.*x + x;
A1 = eps + alpha^3 + (3 * sigma^2 * alpha);
B1 =@(x,lambda) (-3 * lambda * M) * ((hbar(x).^2) + (2 .* hbar(x) .* alpha) + (sigma^2) + (alpha^2));
a1 =@(x) tanh(M .* hbar(x));
b1 =@(x) 1 - ((tanh(M .* hbar(x))).^2);
c1 = (M * alpha) - ((M^3 * A1)/3);
d1 = @(lambda) 2 * (M^2) * (1 + lambda);
C1=@(x) a1(x) + (b1(x) .* c1);
D1 =@(x,lambda) d1(lambda) .* ((hbar(x).^3) + (3 .* (hbar(x).^2) .* alpha) + (3 .* hbar(x) .* (alpha)^2) + (3 .* hbar(x) .* (sigma)^2) + eps + (3 * alpha * (sigma^2)) + (alpha^3));
f1 =@(x,lambda) B1(x,lambda) + (3 * lambda .* C1(x) .* hbar(x)) + (3 * lambda .* C1(x) .* alpha) + (D1(x,lambda) .* C1(x));
f2 =@(x,lambda) 12 * (M^2) * (1 + lambda) .* C1(x);
f3 = psi * (e^3);
f4 = @(lambda) (1 + lambda) *180 * ((1 - e)^2);
f5 =@(lambda) 1/(2 + lambda);
F = @(x,lambda) ((f5(lambda) .* f1(x,lambda))./f2(x,lambda)) + (f3/f4(lambda));
q1 =@(x,lambda) hbar(x) ./ (2 .* F(x,lambda));
Q1 = int(q1,x,0,1); % integration with respect to "x" limit 0 to 1
q2 = @(x,lambda) 1./(F(x,lambda));
Q2 = int(q2,x,0,1); % integration with respect to "x" limit 0 to 1
Q = @(lambda) Q1(lambda)/Q2(lambda); % Q is a function of "lambda"
p1 = @(x,lambda) (1./F(x,lambda)) .* ((0.5 .* hbar(x)) - Q(lambda));
P = @(x,lambda) int(p1,x); % integration with respect to "x" limit 0 to x
% W is obtained by integrating "P" with respect to x limit of x is 0 to
% 1
W =@(lambda) int(P,x,0,1);
fplot(@(lambda) W(lambda))
end
legend(num2str(psi_list'))
ylim([0 1])
xlim([0 4])
set(gca, 'ytick', 0:0.1:1);
set(gca, 'xtick', 0:0.2:4);
xlabel('lambda')
ylabel('W(lambda)')
Sam Chak
on 8 Jul 2022
Can you plot the lines manually one by one? You'll definitely find out what went wrong then.
AVINASH SAHU
on 8 Jul 2022
Accepted Answer
alpha = -0.1;
sigma = 0.1;
eps = -0.1;
e = 0.2;
M = 10;
a = 2;
%figure
psi_list = [10, 20, 30, 40, 50];
lambda = 0:0.2:4;
for i = 1:numel(psi_list)
psi = psi_list(i);
hbar = @(x) a - a.*x + x;
A1 = eps + alpha^3 + (3 * sigma^2 * alpha);
B1 =@(x,lambda) (-3 * lambda * M) * ((hbar(x).^2) + (2 .* hbar(x) .* alpha) + (sigma^2) + (alpha^2));
a1 =@(x) tanh(M .* hbar(x));
b1 =@(x) 1 - ((tanh(M .* hbar(x))).^2);
c1 = (M * alpha) - ((M^3 * A1)/3);
d1 = @(lambda) 2 * (M^2) * (1 + lambda);
C1=@(x) a1(x) + (b1(x) .* c1);
D1 =@(x,lambda) d1(lambda) .* ((hbar(x).^3) + (3 .* (hbar(x).^2) .* alpha) + (3 .* hbar(x) .* (alpha)^2) + (3 .* hbar(x) .* (sigma)^2) + eps + (3 * alpha * (sigma^2)) + (alpha^3));
f1 =@(x,lambda) B1(x,lambda) + (3 * lambda .* C1(x) .* hbar(x)) + (3 * lambda .* C1(x) .* alpha) + (D1(x,lambda) .* C1(x));
f2 =@(x,lambda) 12 * (M^2) * (1 + lambda) .* C1(x);
f3 = psi * (e^3);
f4 = @(lambda) (1 + lambda) *180 * ((1 - e)^2);
f5 =@(lambda) 1./(2 + lambda);
F = @(x,lambda) ((f5(lambda) .* f1(x,lambda))./f2(x,lambda)) + (f3./f4(lambda));
q1 =@(x,lambda) hbar(x) ./ (2 .* F(x,lambda));
Q1 = @(lambda)integral(@(x)q1(x,lambda),0,1,'ArrayValued',true); % integration with respect to "x" limit 0 to 1
q2 = @(x,lambda) 1./(F(x,lambda));
Q2 = @(lambda)integral(@(x)q2(x,lambda),0,1,'ArrayValued',true); % integration with respect to "x" limit 0 to 1
Q = @(lambda) Q1(lambda)./Q2(lambda); % Q is a function of "lambda"
p1 = @(x,lambda) (1./F(x,lambda)) .* ((0.5 .* hbar(x)) - Q(lambda));
P = @(x,lambda) integral(@(y)p1(y,lambda),0,x,'ArrayValued',true); % integration with respect to "x" limit 0 to x
% W is obtained by integrating "P" with respect to x limit to x is 0 to
% 1
W =@(lambda) integral(@(x)P(x,lambda),0,1,'ArrayValued',true);
%fplot(@(lambda) W(lambda))
plot(lambda,W(lambda))
hold on
end
legend(num2str(psi_list'))
ylim([0 1])
xlim([0 4])
set(gca, 'ytick', 0:0.1:1);
set(gca, 'xtick', 0:0.2:4);
xlabel('lambda')
ylabel('W(lambda)')
1 Comment
AVINASH SAHU
on 8 Jul 2022
More Answers (0)
Categories
Find more on Logical in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
