function [valp vecp] = PutereInv(d, s, h, y, maxIter, tol)

  a = s;

  n = length(d);

  for k = 1:maxIter

  daux = d - h;

   z = Thomas(s, daux, s, y);

   %disp(z);

   vecp = z/norm(z);

   %disp(vecp);

   temp = zeros(1, n);

   temp(1) = d(1)*vecp(1) + s(1)*vecp(2);
   temp(n) = s(n-1)*vecp(n-1) + d(n)*vecp(n);

   for i = 2:n-1
    temp(i) = s(i-1)*vecp(i-1) + vecp(i)*d(i) + s(i)*vecp(i+1);
   endfor

   %disp(temp)

   valp = temp*vecp;

   if(norm(temp' - valp*vecp) < tol)
    break;
   endif

   h = valp;

  endfor

endfunction