How to get the derivate using bvp4c

1 Ansicht (letzte 30 Tage)
Syed Mohiuddin
Syed Mohiuddin am 26 Nov. 2022
Kommentiert: Torsten am 27 Nov. 2022
I have a coupled non-linear differential equations u''-b1*(t')*(u')+(1+b1*t)*[G1*F1*t+G2*F1*p-F3*P]=0; t''-b2*(t')^2+B*F6*(u')^2+(b2-b1)*t*B*F6*(u')^2-b2*b1*B*F6*t^2*(u')^2=0; p''- A*p=0 and the boundary conditions are u=0,t=1+m,p=1+n at y=-1 and u=0,t=1,p=1 at y=1.
The code is:
clc;
p=0.01;
Betaf= 207;
Betas = 17;
Beta = 0.5;
kof = 0.613;
kos = 400;
m = 1;
b2 = 0.5;
b1 = 0.5;
G1 = 5;
G2 = 5;
A = 0.5;
Rhof = 997.1;
Rhos = 8933;
P = 0.5;
n=0.5;
B=0.01;
A1 = (1-p).^2.5;
A2 = 1/(1 + 1/Beta);
A3 = (1-p)+p.*((Rhos.*Betas)./(Rhof.*Betaf));
F1 = A2.*A3;
F3 = A1.*A2;
F4 = (kos + 2*kof - 2*p.*(kof - kos))/(kos + 2*kof + p.*(kof - kos));
F5 = (1 + 1/Beta)./A1;
F6 = F5./F4;
dydx=@(x,y)[y(4);
y(5);
y(6);
(b1.*y(4).*y(5)-(1+b1.*y(2)).*(G1.*F1.*y(2)+G2.*F1.*y(3)-F3.*P));
b2.*y(5).^2-B.*F6.*y(4).^2+(b2-b1).*y(2).*B.*F6.*y(4).^2-b2.*b1.*B.*F6.*y(2).^2*y(4).^2;
Alpha.*y(3)];
BC=@(ya,yb)[ya(1);yb(1);ya(2)-(1+m);yb(2)-1.0;ya(3)-(1+n);yb(3)-1.0];
yinit=[0.01;0.01;0.01;0.01;0.01;0.01];
solinit=bvpinit(linspace(-1,1,50),yinit);
S=bvp4c(dydx,BC,solinit)
I like to know how to get the derivative value du/dy at y= -1 and y=1 from the above program
Please help me to complete my code.

Akzeptierte Antwort

Torsten
Torsten am 26 Nov. 2022
I assumed Alpha = A in your code.
clc;
p=0.01;
Betaf= 207;
Betas = 17;
Beta = 0.5;
kof = 0.613;
kos = 400;
m = 1;
b2 = 0.5;
b1 = 0.5;
G1 = 5;
G2 = 5;
A = 0.5;
Rhof = 997.1;
Rhos = 8933;
P = 0.5;
n=0.5;
B=0.01;
A1 = (1-p).^2.5;
A2 = 1/(1 + 1/Beta);
A3 = (1-p)+p.*((Rhos.*Betas)./(Rhof.*Betaf));
F1 = A2.*A3;
F3 = A1.*A2;
F4 = (kos + 2*kof - 2*p.*(kof - kos))/(kos + 2*kof + p.*(kof - kos));
F5 = (1 + 1/Beta)./A1;
F6 = F5./F4;
dydx=@(x,y)[y(4);
y(5);
y(6);
(b1.*y(4).*y(5)-(1+b1.*y(2)).*(G1.*F1.*y(2)+G2.*F1.*y(3)-F3.*P));
b2.*y(5).^2-B.*F6.*y(4).^2+(b2-b1).*y(2).*B.*F6.*y(4).^2-b2.*b1.*B.*F6.*y(2).^2*y(4).^2;
A*y(3)];
BC=@(ya,yb)[ya(1);yb(1);ya(2)-(1+m);yb(2)-1.0;ya(3)-(1+n);yb(3)-1.0];
yinit=[0.01;0.01;0.01;0.01;0.01;0.01];
solinit=bvpinit(linspace(-1,1,50),yinit);
S=bvp4c(dydx,BC,solinit)
S = struct with fields:
solver: 'bvp4c' x: [-1 -0.9388 -0.8776 -0.8163 -0.7551 -0.6939 -0.6327 -0.5714 -0.5102 -0.4490 -0.3878 -0.3469 -0.3265 -0.3061 -0.2449 -0.1837 -0.1224 -0.0612 6.9389e-18 0.0612 0.1020 0.1429 0.1837 0.2449 0.3061 0.3673 0.4286 0.4898 0.5510 0.6122 0.6735 … ] y: [6×37 double] yp: [6×37 double] stats: [1×1 struct]
plot(S.x,S.y(1,:))
S.y(4,1) % u'(-1)
ans = 9.1263
S.y(4,end) % u'(1)
ans = -5.8619
  4 Kommentare
Syed Mohiuddin
Syed Mohiuddin am 27 Nov. 2022
yeh, i know, but i want to make sure the solution is correct. Thank you very much
Torsten
Torsten am 27 Nov. 2022
So you should come to the result that
S.y(5,1) is dt/dy at y = -1
S.y(5,end) is dt/dy at y = 1
S.y(6,1) is dp/dy at y = -1
S.y(6,end) is dp/dy at y = 1

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Produkte


Version

R2016b

Community Treasure Hunt

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

Start Hunting!

Translated by