How to draw Fig. 4

2 visualizaciones (últimos 30 días)
MINATI
MINATI el 30 de Abr. de 2019
Comentada: MINATI el 15 de Mayo de 2019
function main
Pr=1; L=-1;D=1;
R=0.1;Sc=1;
% D=input('D='); %%D=0.5, 1, 1.5
xa=0;xb=6;
solinit=bvpinit(linspace(xa,xb,101),[0 1 0 1 0 1 0]);
sol=bvp4c(@ode,@bc,solinit);
xint=linspace(xa,xb,101);
sxint=deval(sol,xint);
figure(1)
plot(xint,sxint(2,:));
xlabel('\eta');
ylabel('f^\prime');
hold on
function res=bc(ya,yb)
res=[ya(1); ya(2)-L-D*ya(3); ya(4)-1; ya(7)-1; yb(2)-1; yb(4);yb(6)];
end
function dydx=ode(x,y)
dydx=[y(2); y(3); y(2)^2-y(3)*y(1)-1; y(5); -3*Pr*y(1)*y(5)/(3+4*R); y(7); -Sc*y(1)*y(7)];
end
end
Actually I dont know the concept of DUAL solution.

Respuesta aceptada

Stephan
Stephan el 30 de Abr. de 2019
Editada: Stephan el 30 de Abr. de 2019
Hi,
maybe this is useful:
As shown in the link, the results strongly depend on the initial guesses - here is an example:
[sol1, sol2] = main
function [sol1, sol2] = main
Pr=1;
L=-1;
D=[0.5, 1, 1.5];
R=0.1;
Sc=1;
for k = 1:numel(D)
xa=0;xb=6;
solinit1=bvpinit(linspace(xa,xb,101),[0 1 0 1 0 1 0]);
solinit2=bvpinit(linspace(xa,xb,101),[-1 -1 -1 -1 -1 -1 0]);
sol1=bvp4c(@ode,@bc,solinit1);
sol2=bvp4c(@ode,@bc,solinit2);
figure(1)
plot(sol1.x,sol1.yp(1,:),'-r')
xlabel('\eta');
ylabel('f^\prime');
hold on
plot(sol2.x,sol2.yp(1,:),'-b')
end
hold off
function res=bc(ya,yb)
res=[ya(1); ya(2)-L-D(k)*ya(3); ya(4)-1; ya(7)-1; yb(2)-1; yb(4);yb(6)];
end
function dydx=ode(x,y)
dydx=[y(2); y(3); y(2)^2-y(3)*y(1)-1; y(5); -3*Pr*y(1)*y(5)/(3+4*R); y(7); -Sc*y(1)*y(7)];
end
end
results in:
dual.JPG
Best regards
Stephan
  5 comentarios
MINATI
MINATI el 1 de Mayo de 2019
OK Stephan
Thanks
MINATI
MINATI el 15 de Mayo de 2019
Following a Very bad CYCLONE in our State, I didnt follow your code properly as we had no electricity, no internet.
Thats why came late to the forum.
You had considered
f=y(1,:), f'=yp(1,:)
What is f"
Is it yp(2,:), If yes not matching. And
I have taken
f=y(1), f'=y(2); f"= y(3);
Theta= y(4), Theta'=y(5);
phi=y(6);phi'=y(7)

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre MATLAB en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by