function [x,sol]=upwind(i,e,a,b,n,f)
% [x,sol]=upwind(i,e,a,b,n,f)
% solves -e u''+u'=f, u(a)=0, u(b)=1
% with finite differences
% if i=1, centered scheme
% if i=0, upwind
h=(b-a)/n;
A=(diag(2*ones(1,n-1))-diag(ones(1,n-2),-1)-diag(ones(1,n-2),1))/h/h;
x=(a:h:b)';
g=feval(f,x(2:n));
sol=zeros(n+1,1);
if i==1
  A=e*A+(diag(ones(1,n-2),1)-diag(ones(1,n-2),-1))/(2*h);
  sol(n+1)=1;g(n-1)=g(n-1)+e/h/h-1/(2*h);
elseif i==0
  A=e*A+(diag(ones(1,n-1),0)-diag(ones(1,n-2),-1))/h;
  sol(n+1)=1;g(n-1)=g(n-1)+e/h/h;
else
  error('i does not belong to {0,1}')
end 
sol(2:n)=A\g;
