1 view (last 30 days)

function main

D=1; %L=0;

Pr=1;R=0.1;Sc=1;

xa=0;xb=6;

Lv = [-2.5:0.025:0];

p = [];

for i=1:length(Lv)

L = Lv(i);

fODE = @(x,y) [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)];

BCres= @(ya,yb)[ya(1); ya(2)-L-D*ya(3); ya(4)-1; ya(6)-1; yb(2)-1; yb(4);yb(6)];

xint=linspace(xa,xb,101);

solinit=bvpinit(xint,[0 1 0 1 0 1 0]);

sol=bvp4c(fODE,BCres,solinit);

sxint=deval(sol,xint);

%%WE NEED TO PLOT for

S(i,1)=sxint(3,:);

end

figure(1)

plot(Lv,S,'-','Linewidth',1.5);

xlabel('\bf \lambda');

ylabel('\bf C_{f}');

hold on

end

%%While running the code following ERROR occurs:

Subscripted assignment dimension mismatch.

Error in (line 17)

S(i,1)=sxint(3,:);

Walter Roberson
on 19 May 2019

function all_sxint = main

D=1; %L=0;

Pr=1; R=0.1; Sc=1;

xa=0;xb=6;

Lv = [-2.5:0.025:0];

nLv = length(Lv);

all_sxint = cell(nLv, 1);

S = zeros(nLv, 7);

for i=1:nLv

L = Lv(i);

fODE = @(x,y) [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)];

BCres= @(ya,yb)[ya(1); ya(2)-L-D*ya(3); ya(4)-1; ya(6)-1; yb(2)-1; yb(4);yb(6)];

xint=linspace(xa,xb,101);

solinit=bvpinit(xint,[0 1 0 1 0 1 0]);

sol=bvp4c(fODE,BCres,solinit);

sxint=deval(sol,xint);

all_sxint{i} = sxint;

S(i,:) = sxint(3,:);

end

figure(1)

plot(Lv, S, '-', 'Linewidth', 1.5);

xlabel('\bf \lambda');

ylabel('\bf C_{f}');

legend({'bc1', 'bc2', 'bc3', 'bc4', 'bc5', 'bc6', 'bc7'})

end

Assign the output to a variable so that you can examine all of the time points for all of the Lv values afterwards, as you indicate that you need to be able to do that.

Walter Roberson
on 19 May 2019

I noticed some oddities in the output. I used options to push the permitted mesh points way up, and zoomed in closer. It turns out there is something unusual going on at -2.4648 .

Warning: Unable to meet the tolerance without using more than 500000 mesh points.

The last mesh of 292855 points and the solution are available in the output argument.

The maximum residual is 0.204915, while requested accuracy is 0.001.

The above was for Lv = linspace(-2.475, -2.425, 50);

The residual is staying pretty much the same as I increase the number of mesh points and zoom in more closely, which is potentially hinting at a singularity.

Different lines correspond to different xint.

Matt J
on 18 May 2019

Edited: Matt J
on 18 May 2019

sxint(3,:) is not a scalar, but the left hand side S(i,1) is a scalar location.

Walter Roberson
on 19 May 2019

Opportunities for recent engineering grads.

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

Start Hunting!
## 0 Comments

Sign in to comment.