Wrong graph is coming with this code

Question

MINATI 2020-8-7

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/576634-wrong-graph-is-coming-with-this-code

评论： Star Strider 2020-8-27

Graph.pdf

function main
A=0.5; pr=1; a=1; phi=0.1;rhos=997.1;Cps=4179;ks=0.613;rhof=8933;Cpf=385;kf=401; 
a1=(1-phi)^-2.5*((1-phi)+phi*(rhos/rhof));  Knf=(ks+2*kf-2*phi*(kf-ks))/(ks+2*kf+phi*(kf-ks));
a2=((1-phi)+phi*((rhos*Cps)/(rhof*Cpf)))*Knf;
xa = 0;xb = 5;solinit=bvpinit(linspace(xa,xb,101),[0 1 0 a 1 0 0 0]);
sol = bvp4c(@ode,@bc,solinit); x = linspace(xa,xb,101);S = deval(sol,x);
function res = bc(ya,yb)
res = [ya(1); ya(2)-1; ya(4); ya(5)-a; ya(7)-1; yb(2); yb(5); yb(7)];
end 
function dydx = ode(~,y)    
dydx = [y(2);  y(3);  2*a1*y(2)*(y(2)+y(5))-a1*y(3)*(y(1)+y(4));  
        y(5);  y(6);  2*a1*y(5)*(y(2)+y(5))-a1*y(6)*(y(1)+y(4));
        y(8);  A*pr*a2*y(7)*(y(2)+y(5))-pr*a2*y(8)*(y(1)+y(4))];              
end
figure(14)
f1 = @(x,y) S(2);
fsurf(f1)
xlabel '\bfx';ylabel '\bfy';zlabel '\bff^\prime(\eta)'
hold on
end
%% The graph is not obeying Boundary Condition
%% I am trying graphs like Fig. (14)  and  another  is contour (Attached)

2 个评论
显示无隐藏无

MINATI 2020-8-10

can anyone share idea here

MINATI 2020-8-26

Dear Star Strider

Please have a try

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2020-8-26

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/576634-wrong-graph-is-coming-with-this-code#answer_485315

在 MATLAB Online 中打开

Using fsurf is not appropriate here.

Try this instead:

figure(14)
hs = surf(S);
grid on
xlabel '\bfx';ylabel '\bfy';zlabel '\bff^\prime(\eta)'
% shading('interp')                                               % Optional
% hs.EdgeColor = 'none';                                          % Optional

That produces this plot:

I am not certain what the reference to ‘S(2)’ is suposed to do, beccause ‘S’ is a matrix. This plots the entire matrix.

.

6 个评论
显示 4更早的评论隐藏 4更早的评论

Star Strider 2020-8-26

在 MATLAB Online 中打开

My pleasure!

I am completely lost.

With respect to: ‘How to draw surf(S(2,:)) with x & y ’ this works:

hs = surf(x, [2 2], [1;1]*S(2,:))

With respect to: ‘And also surf(meshgrid(x), meshgrid(phi), S(3,1))’ is just not going to work because ‘S(3,1)’ is a point, so plotting it is useless, as welll as the rest of that not making any sense.

With respect to: ‘Where phi = linspace(0,0.15,101)’ note that ‘phi’ is currently a scalar. I have no idea what you are doing, however you will need to incorporate that vector with your code. Every variable that uses ‘phi’ will also become vectors, likely not something bvp4c can handle, since your ODE/BVP function needs to return a column vector.

So the best way of dealing with that may be a loop, looping through the various values of ‘phi’ and storing the results by concatenating the resulting ‘S’ matrices along the third dimension. Then choose the appropriate values (rows, columns, pages) of that 3D matrix for the plots.

When you have that working, post back here with the results. It will likely be possible to plot them.

Star Strider 2020-8-26

MINATI’s Answer moved here —

Dear Star Strider have a look into the following code:

function main

A = 0.5; pr = 1; %%a = 1; phi = 0.1;

rhos = 997.1;Cps=4179;ks=0.613;rhof=8933;Cpf=385;kf=401;

phiv = linspace (0,0.15,101);

for a = [0 0.2 0.5 1]

for i = 1:length(phiv)

phi = phiv(i);

a1=(1-phi)^2.5*((1-phi)+phi*(rhos/rhof)); Knf=(ks+2*kf-2*phi*(kf-ks))/(ks+2*kf+phi*(kf-ks));

a2 =((1-phi)+phi*((rhos*Cps)/(rhof*Cpf)))*Knf;

BCres = @(ya,yb)[ya(1); ya(2)-1; ya(4); ya(5)-a; ya(7)-1; yb(2); yb(5); yb(7)];

fODE = @(x,y)[y(2); y(3); 2*a1*y(2)*(y(2)+y(5))-a1*y(3)*(y(1)+y(4));

y(5); y(6); 2*a1*y(5)*(y(2)+y(5))-a1*y(6)*(y(1)+y(4));

y(8); A*pr*a2*y(7)*(y(2)+y(5))-pr*a2*y(8)*(y(1)+y(4))];

xa = 0;xb = 5;solinit=bvpinit(linspace(xa,xb,101),[0 1 0 a 1 0 0 0]);

sol = bvp4c(fODE,BCres,solinit);x = linspace(xa,xb,101);S = deval(sol,x);

Sk(i,1) = S(3,1)/((1-phi)^2.5); Nu(i,1) = - Knf * S(5,1);

end

figure(10)

plot(phiv,Sk,'-','Linewidth',1.5)

% % title '(b)'

xlabel('\bf\phi');

ylabel('\bfSk_{r}')

hold on

end

%%%%%%%%I have incorporated a loop for phi variation where plotting of S(3,1) is OK but

surf(S(3,:)) is not coming

Can it possible to modify for surf of Sk (OR In x-axis, phi variation In Y-axis for a = linspace(0,1,101))

Please do something

Star Strider 2020-8-26

在 MATLAB Online 中打开

Create a second loop for ‘a’:

function main
A = 0.5; pr = 1; %%a = 1; phi = 0.1;
rhos = 997.1;Cps=4179;ks=0.613;rhof=8933;Cpf=385;kf=401; 
phiv = linspace (0,0.15,101);
av = [0 0.2 0.5 1];
for k = 1:numel(av)
    a = av(k);
 for i = 1:length(phiv)
     phi = phiv(i);
a1=(1-phi)^2.5*((1-phi)+phi*(rhos/rhof));  Knf=(ks+2*kf-2*phi*(kf-ks))/(ks+2*kf+phi*(kf-ks));
a2 =((1-phi)+phi*((rhos*Cps)/(rhof*Cpf)))*Knf;
BCres = @(ya,yb)[ya(1); ya(2)-1; ya(4); ya(5)-a; ya(7)-1; yb(2); yb(5); yb(7)];  
fODE = @(x,y)[y(2);  y(3);  2*a1*y(2)*(y(2)+y(5))-a1*y(3)*(y(1)+y(4));  
        y(5);  y(6);  2*a1*y(5)*(y(2)+y(5))-a1*y(6)*(y(1)+y(4));
        y(8);  A*pr*a2*y(7)*(y(2)+y(5))-pr*a2*y(8)*(y(1)+y(4))]; 
xa = 0;xb = 5;solinit=bvpinit(linspace(xa,xb,101),[0 1 0 a 1 0 0 0]);
sol = bvp4c(fODE,BCres,solinit);x = linspace(xa,xb,101);S = deval(sol,x);
Sk(i,k) = S(3,1)/((1-phi)^2.5);   Nu(i,1) =  - Knf * S(5,1);  
 end
end
figure(10)
surf(av, phiv, Sk)
% % title '(b)'
xlabel('a')
ylabel('\bf\phi'); 
zlabel('\bfSk_{r}')
view(40,40)
end

Also, put the surf call after both loops.

This then creates a matrix for ‘Sk’ that you can use with surf.

.

MINATI 2020-8-27

Thanks Star Strider

You are awsome.

Any way the first problem remain unanswered (Attached pdf)

surf(x, y, S(2,:)) and contour

pdf exist means code can be created but I can't.

But I have confidence on You in this matter. If you agree take your time to create.

A great thanks!

Star Strider 2020-8-27

Thank you!

Your ccode does not produce the plot in the .PDF file. I have no idea what you are doing.

You can use essentially the same code to create the contour plot as the surf plot. Just call contour instead of surf, and remove the zlabel and view calls.

As always, my pleasure!

请先登录，再进行评论。

Wrong graph is coming with this code

2 个评论
显示无隐藏无

采纳的回答

6 个评论
显示 4更早的评论隐藏 4更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

Wrong graph is coming with this code

2 个评论 显示 无隐藏 无

采纳的回答

6 个评论 显示 4更早的评论隐藏 4更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

2 个评论
显示无隐藏无

6 个评论
显示 4更早的评论隐藏 4更早的评论