How can I solve a system of ODEs having coefficients in vector form using bvp4c ?

조회 수: 5(최근 30일)
Tanya Sharma 3 Oct 2019
댓글: darova 11 Oct 2019
Hi,
I am solving a system of 5 ODEs that contains known parameters, unknown parameters and some coefficients are of the form of vector.
%-------------------------------ODE system-----------------------------------%
function eqns = odes(x,y,e)
global Pr phi Ra Da Fr A1 A2 fdesh fdeshdesh thetadesh;
eqns = [y(2)
y(3)
(phi./Da).*y(2)+(2.*phi.*Fr./A1).*fdesh.*y(2)-(fdesh.*1./A1).*y(3)-(fdeshdesh.*1./A1).*y(1)+(2.*fdesh.*1./A1).*y(2)-(e./A1).*y(2)-(phi.*Ra./(A1^2).*A2).*y(4)
y(5)
end
%------------------------------------------------------------------------------------------------
--------------------% fdesh, fdeshdesh and thetadesh are of vector form. these are the known solution of the governing equations .%----------------------
%----------------------------------------------------------------------------------------------
%--------------------------boundary conditions-----------------------------%
function res = ode_bc(ya,yb,e)
res = [ya(1)
ya(2)
ya(3)
ya(4)
yb(2)
yb(4)];
end
%-------------------------------------------------------------------------------------------------
I am getting this error:
Error using bvparguments (line 108)
Error in calling BVP4C(ODEFUN,BCFUN,SOLINIT):
The derivative function ODEFUN should return a column vector of length 5.
%-------------------------------------------------------------------------------------------------------
are the known solutions causing a problem to solve the system?

댓글을 달려면 로그인하십시오.

채택된 답변

darova 3 Oct 2019
Try this
function eqns = odes(x,y,e,x0)
global Pr phi Ra Da Fr A1 A2 fdesh fdeshdesh thetadesh; % global variables are not recommended
% x0 - vector of corresponding values for fdesh, fdeshdesh, thetadesh
% x0(end) should not be bigger than x(end)
fd = interp1(x0,fdesh,x);
fdd = interp1(x0,fdeshdesh,x);
eqns = [y(2)
y(3)
(phi./Da).*y(2)+(2.*phi.*Fr./A1).*fd.*y(2)-(fd.*1./A1).*y(3)-(fdd.*1./A1).*y(1)+(2.*fd.*1./A1).*y(2)-(e./A1).*y(2)-(phi.*Ra./(A1^2).*A2).*y(4)
y(5)
-(Pr./A2).*(fd.*y(5)+thd.*y(1)+e.*y(4))];
end
댓글 수: 10표시숨기기 이전 댓글 수: 9
darova 9 Oct 2019
I don't think β patameter can be found. In the paper you attached everywhere is said that it can be obtained with guess

댓글을 달려면 로그인하십시오.

추가 답변(1개)

Tanya Sharma 10 Oct 2019
Yes and the guess is provided in the "solinit" structure as an extra parameter
beta = 99;
init2 = bvpinit(linspace(-1,1,50),@math4init,beta);
Which we can later be obtain by "sol.parameters", this provides us the value bvp4c has taken nearer to our provided guess.
Which in this case is:
beta=77.8263.
댓글 수: 3표시숨기기 이전 댓글 수: 2
darova 11 Oct 2019
Cool. It works!

댓글을 달려면 로그인하십시오.

R2018b

Community Treasure Hunt

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

Start Hunting!

Translated by