function bvppractice
clear all

xlow=0;
xhigh=1;
solinit = bvpinit(linspace(xlow,xhigh,20),[1 0]);
sol = bvp4c(@bvp4ode,@bvp4bc,solinit);
xint = linspace(xlow,xhigh);
Sxint = deval(sol,xint);
ep=0.1;
analyt=(exp(xint/ep)-1)/(exp(1/ep)-1);
plot(xint,Sxint(1,:),xint,analyt,'r')
xlabel('x')
ylabel('y')
errorvec=abs((analyt-Sxint(1,:))./analyt)
figure
plot(xint,errorvec)
ylabel('error')
xlabel('x')
disp('press any key to continue')
pause

% chip cooling fins
xlow=0;
xhigh=0.025;
solinit = bvpinit(linspace(xlow,xhigh,20),[30 0]);
sol = bvp4c(@bvppin,@bvppinbc,solinit);
xint = linspace(xlow,xhigh);
Sxint = deval(sol,xint);
plot(xint,Sxint(1,:))
xlabel('x (m)')
ylabel('Temperature (C)')
disp('press any key to continue')
pause

% resistor temperature
xlow=1e-8;
xhigh=0.001;
solinit = bvpinit(linspace(xlow,xhigh,20),[30 0]);
sol = bvp4c(@bvpres,@bvpresbc,solinit);
xint = linspace(xlow,xhigh);
Sxint = deval(sol,xint);
plot(xint,Sxint(1,:))
xlabel('r (m)')
ylabel('Temperature (C)')
disp('press any key to continue')
pause

% resistor temperature with convection
xlow=1e-8;
xhigh=0.001;
solinit = bvpinit(linspace(xlow,xhigh,20),[30 0]);
sol = bvp4c(@bvpres2,@bvpres2bc,solinit);
xint = linspace(xlow,xhigh);
Sxint = deval(sol,xint);
plot(xint,Sxint(1,:))
xlabel('r (m)')
ylabel('Temperature (C)')

% -----------------------------------------------
function dydx = bvp4ode(x,y)
eps=0.1;
dydx = [  y(2)
         y(2)/eps ];
% -----------------------------------------------
function res = bvp4bc(ya,yb)
res = [  ya(1)
         yb(1)-1 ];
% -----------------------------------------------
function dydx = bvppin(x,y)
hcka=4000;
tf=20;
dydx = [  y(2)
         hcka*(y(1)-tf) ];
% -----------------------------------------------
function res = bvppinbc(ya,yb)
res = [  ya(1)-40
         yb(2) ];
%
function dydx = bvpres(r,y)
k=0.1;
q=2.1e6;
dydx = [  y(2)
         -q/k-y(2)/r ];
% -----------------------------------------------
function res = bvpresbc(ya,yb)
res = [  ya(2)
         yb(1)-20 ];
%
function dydx = bvpres2(r,y)
k=0.1;
q=2.1e6;
dydx = [  y(2)
         -q/k-y(2)/r ];
% -----------------------------------------------
function res = bvpres2bc(ya,yb)
h=10;
Te=20;
R=0.001;
q=2.1e6;
res = [  ya(2)
         h*(yb(1)-Te)-q*R/2 ];