eps=input('  eps=  ');

    betaL=1; betaR=2;

for NN=1:6

    N=2^(NN+4)

    h=2/N; 

    j=[1:1:N-1];  x=-1+j*h;

    for i=1:N-1
       s=x(i); p1=-s;  q1=-1;
       s1=(s+1)/eps;  s2=(s-1)/eps;
       u_true(i)=exp(-s1)+2*exp(s2);
       f1=(s1-1)*exp(-s1)-2*(s2+1)*exp(s2);

          for j=1:N-1
            if i==j
              A(i,j)=-2*eps+q1*h^2;
            elseif i==j+1
              A(i,j)=eps-0.5*h*p1;
            elseif j==i+1
              A(i,j)=eps+0.5*h*p1;
            else
              A(i,j)=0;
            end
          end
      b(i)=f1*h^2;
    end

       ss1=eps+0.5*h*x(N-1);
       b(N-1)=b(N-1)-ss1*betaR;

       ss2=eps-0.5*h*x(1);
       b(1)=b(1)-ss2*betaL;

    u=A\b';

    for i=1:N-1
       error(i)=abs(u(i)-u_true(i));
    end

    xx(NN)=N;
    err(NN)=max(error);

    end

    fprintf(1,      '%16.0f          %13.3e  \n', [xx;  err]);