Filter löschen
Filter löschen

how to write matlab bvp4c code for moving boundary surface problem?

9 Ansichten (letzte 30 Tage)
hello, I am new in the bvp4c method of Matlab. I want to find multiple graphs for different value of 'r' for the following problem. Equation is f''' +ff''=0 Boundary condition is f(0)=0, f'(0)=1-r, f'(infnty)=r

Akzeptierte Antwort

Stephan
Stephan am 6 Nov. 2018
Bearbeitet: Stephan am 6 Nov. 2018
Hi,
i recommend to start simple and work step by step:
  1. Have a read in the documentation of bvp4c
  2. Try to implement your model with a single fixed value for r and get a valid solution
  3. Once you have this, use a for loop to vary r like needed for your purposes and save all the results in a Matrix in different columns or rows
  4. Plot the results by using the saved values
Problems by doing this? --> Come back and ask a specific question including your code so far.
Best regards
Stephan
  2 Kommentare
Anil Gautam
Anil Gautam am 7 Nov. 2018
Bearbeitet: Torsten am 7 Nov. 2018
function cortell
global r
r=[0]
solinit = bvpinit(linspace(0,10,1000),[0 1- r r]);
options = bvpset('stats','on');
sol = bvp4c(@equation,@bvpbc,solinit,options);
eta = sol.x;
u = sol.y;
clf reset
plot(eta,u(3,:));
hold on;
axis([0 4 0 1.3]);
xlabel('\eta')
ylabel('f''(\eta)')
shg
% --------------------------------------------------------------------------
function dy=equation(~,y)
%dy=zeros(4,1);
dy=[y(2);
y(3);
-y(1)*y(3);
%---------------------------------------------------------------------
];
function res = bvpbc(y0,yinf)
global r
res = [y0(1)
y0(2)-1+ r
yinf(2)- r
];
if true
% code
end
this is code . how I can I get the different graph for different value of 'r' . In other word, for r=0.0, 0.25,0.75 how i can plot graph for f on y axis and eta on x axis.. please suggest.
Stephan
Stephan am 7 Nov. 2018
See the answer from Torsten...

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (1)

Torsten
Torsten am 7 Nov. 2018
Bearbeitet: Torsten am 7 Nov. 2018
function main
r=[0 0.02 0.04]
for i=1:numel(r)
r_actual = r(i);
solinit = bvpinit(linspace(0,10,10000),[0 1-r_actual 0]);
options = bvpset('stats','on');
sol{i} = bvp4c(@equation,@(y0,yinf)bvpbc(y0,yinf,r_actual),solinit,options);
end
for i=1:numel(r)
eta = sol{i}.x;
u = sol{i}.y;
plot(eta,u(3,:));
hold on;
end
end
% --------------------------------------------------------------------------
function dy=equation(~,y)
%dy=zeros(4,1);
dy=[y(2);
y(3);
-y(1)*y(3);
%---------------------------------------------------------------------
];
end
function res = bvpbc(y0,yinf,r_actual)
res = [y0(1)
y0(2)-1+ r_actual
yinf(2)- r_actual
];
end
  3 Kommentare
Torsten
Torsten am 7 Nov. 2018
Bearbeitet: Torsten am 7 Nov. 2018
Try this code - it uses the solution of the last step as initial guess for the next:
function main
global ix
r=[0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0];
r_actual = r(1);
ir=1;
solinit = bvpinit(linspace(0,10,1000),@(x)guess(x,ir,[0 1-r_actual r_actual]));
options = bvpset('stats','on');
ix = 0;
sol{1} = bvp4c(@equation,@(y0,yinf)bvpbc(y0,yinf,r_actual),solinit,options);
for ir=2:numel(r)
r_actual = r(ir);
solinit = bvpinit(sol{ir-1}.x,@(x)guess(x,ir,sol{ir-1}.y));
options = bvpset('stats','on');
ix = 0;
sol{ir} = bvp4c(@equation,@(y0,yinf)bvpbc(y0,yinf,r_actual),solinit,options);
end
for ir=1:numel(r)
eta = sol{ir}.x;
u = sol{ir}.y;
plot(eta,u(3,:));
hold on;
end
end
% --------------------------------------------------------------------------
function dy=equation(~,y)
%dy=zeros(4,1);
dy=[y(2);
y(3);
-y(1)*y(3);
%---------------------------------------------------------------------
];
end
function res = bvpbc(y0,yinf,r_actual)
res = [y0(1)
y0(2)-1+ r_actual
yinf(2)- r_actual
];
end
function yini=guess(x,ir,sol)
global ix
ix = ix + 1;
if ir==1
yini(1)=sol(1);
yini(2)=sol(2);
yini(3)=sol(3);
else
yini(1)=sol(1,ix);
yini(2)=sol(2,ix);
yini(3)=sol(3,ix);
end
end

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Specifying Target for Graphics Output finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by