Plot a 2d figure for a complicated function

1 view (last 30 days)
Hi. I wanted to plot a 2d figure for the equation below (lambda) that are a function of c, but it doesn't give me any figures. Could you please tell me what is my problem? Thanks in advance for any help
vplc=0.25;delta=2.5;tau_max=44000;Ktau=0.045;%tauP=0.027;
kc=0.1; kh=0.05;Vp=0.9;Kbar=0.000015;kp=0.15;gamma=5.5;kb=0.4;
Vpm=0.000159;Kpm=0.15;alpha0=.00000681;alpha1=0.0000227;Ke=7;vs=0.002;ks=0.1;
Kf=0.18;kplc=0.055;ki=2;gamma=5.5;kipr=0.18;
[c]=meshgrid(0.001:0.005:1);
%[c]=meshgrid(0.002:0.005:1);
%c=0.07362169;h=0.5823305;ce=16.64219;p=0;
p=(vplc./ki).*(c.^2./((kplc).^2+c.^2));
h=kh.^4./(c.^4+kh.^4);
Po=(p.^2.*c.^4.*h)./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4));
A=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^4.*c.^8.*h.*(1+kb))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
B=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^2.*c.^4.*h.*(4.*c.^3.*p.^2.*h.*(1+kb)+4.*c.^3.*kb.*kp.^2))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
T=kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)-c.^4.*(1./tau_max);
D=(kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)).*(-c.^4.*(1./tau_max))-((kipr.*B.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))).*(c.^4.*(1./tau_max).*((-4.*c.^3.*kh.^4)./(c.^4+kh.^4).^2)));
lambda=(T+sqrt(T.^2-4.*D))./2;
plot(c,lambda)
hold on

Accepted Answer

Voss
Voss on 20 Jul 2022
vplc=0.25;delta=2.5;tau_max=44000;Ktau=0.045;%tauP=0.027;
kc=0.1; kh=0.05;Vp=0.9;Kbar=0.000015;kp=0.15;gamma=5.5;kb=0.4;
Vpm=0.000159;Kpm=0.15;alpha0=.00000681;alpha1=0.0000227;Ke=7;vs=0.002;ks=0.1;
Kf=0.18;kplc=0.055;ki=2;gamma=5.5;kipr=0.18;
c is a matrix where every row is identical
[c]=meshgrid(0.001:0.005:1)
c = 200×200
0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460 0.0010 0.0060 0.0110 0.0160 0.0210 0.0260 0.0310 0.0360 0.0410 0.0460 0.0510 0.0560 0.0610 0.0660 0.0710 0.0760 0.0810 0.0860 0.0910 0.0960 0.1010 0.1060 0.1110 0.1160 0.1210 0.1260 0.1310 0.1360 0.1410 0.1460
% is c equivalent to repeated copies of its first row?
isequal(c,repmat(c(1,:),size(c,1),1))
ans = logical
1
% c = 0.001:0.005:1;
%[c]=meshgrid(0.002:0.005:1);
%c=0.07362169;h=0.5823305;ce=16.64219;p=0;
p=(vplc./ki).*(c.^2./((kplc).^2+c.^2));
h=kh.^4./(c.^4+kh.^4);
Po=(p.^2.*c.^4.*h)./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4));
A=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^4.*c.^8.*h.*(1+kb))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
B=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^2.*c.^4.*h.*(4.*c.^3.*p.^2.*h.*(1+kb)+4.*c.^3.*kb.*kp.^2))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
T=kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)-c.^4.*(1./tau_max);
D=(kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)).*(-c.^4.*(1./tau_max))-((kipr.*B.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))).*(c.^4.*(1./tau_max).*((-4.*c.^3.*kh.^4)./(c.^4+kh.^4).^2)));
lambda=(T+sqrt(T.^2-4.*D))./2;
so lambda is also a matrix with every row identical
% is lambda equivalent to repeated copies of its first row?
isequal(lambda,repmat(lambda(1,:),size(lambda,1),1))
ans = logical
1
Using plot with matrices, one line is created for each column of the matrices
h = plot(c,lambda) % 200 lines plotted
h =
200×1 Line array: Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line
So when every row is the same, you get lines consisting of multiple copies (200 in this case) of the same point.
get(h(1),'XData') % x-coordinates of the points in the first line
ans = 1×200
1.0e-03 * 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000
get(h(1),'YData') % y-coordinates of the points in the first line
ans = 1×200
1.0e-16 * -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274 -0.2274
You can use a data marker to see the points/lines:
h = plot(c,lambda,'.') % different colored points for different lines
h =
200×1 Line array: Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line Line
Therefore, since c is redundant, in that it is a matrix consisting of 200 copies of the same row, perhaps it should be a vector?
c = 0.001:0.005:1;
%c=0.07362169;h=0.5823305;ce=16.64219;p=0;
p=(vplc./ki).*(c.^2./((kplc).^2+c.^2));
h=kh.^4./(c.^4+kh.^4);
Po=(p.^2.*c.^4.*h)./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4));
A=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^4.*c.^8.*h.*(1+kb))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
B=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^2.*c.^4.*h.*(4.*c.^3.*p.^2.*h.*(1+kb)+4.*c.^3.*kb.*kp.^2))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
T=kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)-c.^4.*(1./tau_max);
D=(kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)).*(-c.^4.*(1./tau_max))-((kipr.*B.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))).*(c.^4.*(1./tau_max).*((-4.*c.^3.*kh.^4)./(c.^4+kh.^4).^2)));
lambda=(T+sqrt(T.^2-4.*D))./2;
In which case lambda, T and D are also vectors of the same size, so plotting them is no problem:
figure % using different widths and styles to distinguish the lines:
plot(c,lambda,'LineWidth',3) % thick blue line
hold on
plot(c,T) % red line
hold on
plot(c,D,'--','LineWidth',2) % dashed yellow line

More Answers (1)

Chunru
Chunru on 20 Jul 2022
vplc=0.25;delta=2.5;tau_max=44000;Ktau=0.045;%tauP=0.027;
kc=0.1; kh=0.05;Vp=0.9;Kbar=0.000015;kp=0.15;gamma=5.5;kb=0.4;
Vpm=0.000159;Kpm=0.15;alpha0=.00000681;alpha1=0.0000227;Ke=7;vs=0.002;ks=0.1;
Kf=0.18;kplc=0.055;ki=2;gamma=5.5;kipr=0.18;
% grid
x = 0.001:0.005:1;
y = 0.001:0.005:1;
[c]=meshgrid(x, y);
%[c]=meshgrid(0.002:0.005:1);
%c=0.07362169;h=0.5823305;ce=16.64219;p=0;
p=(vplc./ki).*(c.^2./((kplc).^2+c.^2));
h=kh.^4./(c.^4+kh.^4);
Po=(p.^2.*c.^4.*h)./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4));
A=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^4.*c.^8.*h.*(1+kb))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
B=(4.*c.^3.*p.^2.*h.*(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^4+c.^4))-p.^2.*c.^4.*h.*(4.*c.^3.*p.^2.*h.*(1+kb)+4.*c.^3.*kb.*kp.^2))./(p.^2.*(1+kb).*c.^4.*h+kb.*kp.^2.*(kc.^2+c.^4)).^2;
T=kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)-c.^4.*(1./tau_max);
D=(kipr.*A.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))-kipr.*Po.*(1+gamma)-((2.*c.*vs.*(c.^2+ks.^2)-2.*vs.*c.^3)./(c.^2+ks.^2).^2)).*(-c.^4.*(1./tau_max))-((kipr.*B.*((vs.*c.^2./(c.^2+ks.^2)).*(1./kipr.*Po))).*(c.^4.*(1./tau_max).*((-4.*c.^3.*kh.^4)./(c.^4+kh.^4).^2)));
lambda=(T+sqrt(T.^2-4.*D))./2;
%whos
% This is a 3d plot rather than 2d plot
%plot(c,lambda)
surf(x, y, lambda, 'EdgeColor', 'none')
hold on

Community Treasure Hunt

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

Start Hunting!