stoppingCriteria <- function(x, y, TOL = .Machine$double.eps, method = c("1", "2", "3"), theta = NULL) {
  method <- match.arg(method)

  if (method == "1") {
    return(abs(y - x) < TOL)
  } else if (method == "2") {
    if (y == 0)
      stop("This choice of stopping criteria not good for solutions near zero")
    return((abs(y - x) / abs(y)) < TOL)
  } else {
    if (is.null(theta))
      stop("theta > 0 needs to be specified if using method 3")
    return(abs(y - x) / max(abs(y), theta) < TOL)
  }
}