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)
plot(S.x,S.y(1,:))
S.y(4,1) % u'(-1)
S.y(4,end) % u'(1)
