MATLAB: How to find a solution for this non-linear, stiff 2nd order ODE

nonlinearnumerical integrationode

I need to find a solution to an equation for T in the form:
aT'' – bT^4 + c = 0, T(0) = T0, T'(L) = 0
(FYI: this is a steady-state heat conduction and emission problem along a cantilever. The clamped end of the rod is kept at T(0)=T0, and there is no conduction at the loose end T'(L)=0. In a vacuum, so heat is only dissipated via black-body emission)
My code:
% Physical parameters (fixed):
a = 0.0018;
b = 4.2525e-08;
c = 11.3264;
L = 1.5;
T0 = 300;
rinit = 100; % initial guess parameter
odefun = @(x,y) [ y(2); (a/b)*y(1)^4 - c/a ];
bcfun = @(ya,yb) [ya(1)-T0; yb(2)];
xinit = linspace(0,L,1e3);
yinit = @(x) [rinit*x.^2 - 2*rinit*L*x + T0; 2*rinit*x - 2*rinit*L]; % quadratic initial guess for the solution
plot(xinit, yinit(xinit))
solinit=bvpinit(xinit, yinit);
sol=bvp4c(odefun,bcfun,solinit);
plot(sol.x,sol.y(1,:))
If you run it like that, it does not compute a solution. If you run it with
a = 1;
b = 1;
c = 0;
L = 1;
T0 = 1;
rinit = 1;
you should see what kind of solution I am looking for.

Best Answer

odefun = @(x,y) [ y(2); b/a*y(1)^4 - c/a ];
Best wishes
Torsten.