#include <bits/stdc++.h>
using namespace std;
template <typename A, typename B>
inline void chmin(A &a, B b) {
  if (a > b) a = b;
}
template <typename A, typename B>
inline void chmax(A &a, B b) {
  if (a < b) a = b;
}
template <int64_t mod>
struct ModInt {
  int64_t a;
  inline int64_t normalize(int64_t x) {
    x %= mod;
    if (x < 0) x += mod;
    return x;
  }
  ModInt(const int64_t a = 0) : a(normalize(a)) {}
  ModInt &operator+=(const ModInt &x) {
    a += x.a;
    if (a >= mod) a -= mod;
    return *this;
  }
  ModInt &operator-=(const ModInt &x) {
    a -= x.a;
    if (a < 0) a += mod;
    return *this;
  }
  ModInt &operator*=(const ModInt &x) {
    a = a * x.a % mod;
    return *this;
  }
  ModInt &operator/=(const ModInt &x) {
    *this *= x.inv();
    return *this;
  }
  ModInt operator+(const ModInt x) { return ModInt(*this) += x; }
  ModInt operator-(const ModInt x) { return ModInt(*this) -= x; }
  ModInt operator*(const ModInt x) { return ModInt(*this) *= x; }
  ModInt operator/(const ModInt x) { return ModInt(*this) /= x; }
  bool operator==(const ModInt x) { return a == x.a; }
  bool operator!=(const ModInt x) { return a != x.a; }
  ModInt operator-() { return ModInt(0) - ModInt(*this); }
  ModInt pow(int64_t ex) const {
    int64_t x = a;
    int64_t res = 1;
    while (ex) {
      if (ex & 1) res = res * x % mod;
      x = x * x % mod;
      ex >>= 1;
    }
    return ModInt(res);
  }
  ModInt inv() const { return pow(mod - 2); }
};
template <int64_t mod>
istream &operator>>(istream &in, ModInt<mod> &a) {
  return (in >> a.a);
}
template <int64_t mod>
ostream &operator<<(ostream &out, const ModInt<mod> &a) {
  return (out << a.a);
}
using mint = ModInt<998244353>;
template <class T>
struct Polynomial : vector<T> {
  using vector<T>::vector;
  Polynomial &operator+=(const Polynomial &x) {
    if (vector<T>::size() < x.size()) vector<T>::resize(x.size());
    for (long long i = 0; i < x.size(); i++) vector<T>::operator[](i) += x[i];
    return *this;
  }
  Polynomial &operator-=(const Polynomial &x) {
    if (vector<T>::size() < x.size()) vector<T>::resize(x.size());
    for (long long i = 0; i < x.size(); i++) vector<T>::operator[](i) -= x[i];
    return *this;
  }
  Polynomial &operator*=(const Polynomial &x) {
    if (vector<T>::size() == 0 || x.size() == 0) {
      vector<T>::clear();
      return *this;
    }
    Polynomial z(vector<T>::size() + x.size() - 1);
    for (long long i = 0; i < vector<T>::size(); i++) {
      for (long long j = 0; j < x.size(); j++) {
        z[i + j] += vector<T>::operator[](i) * x[j];
      }
    }
    return *this = z;
  }
  Polynomial &operator+=(const T &x) {
    if (vector<T>::size() == 0) vector<T>::resize(1);
    vector<T>::operator[](0) += x;
    return *this;
  }
  Polynomial &operator-=(const T &x) {
    if (vector<T>::size() == 0) vector<T>::resize(1);
    vector<T>::operator[](0) -= x;
    return *this;
  }
  Polynomial &operator*=(const T &x) {
    for (long long i = 0; i < vector<T>::size(); i++)
      vector<T>::operator[](i) *= x;
    return *this;
  }
  Polynomial &operator/=(const T &x) {
    for (long long i = 0; i < vector<T>::size(); i++)
      vector<T>::operator[](i) /= x;
    return *this;
  }
  Polynomial operator+(const Polynomial &x) const {
    return Polynomial(*this) += x;
  }
  Polynomial operator-(const Polynomial &x) const {
    return Polynomial(*this) -= x;
  }
  Polynomial operator*(const Polynomial &x) const {
    return Polynomial(*this) *= x;
  }
  Polynomial operator+(const T &x) const { return Polynomial(*this) += x; }
  Polynomial operator-(const T &x) const { return Polynomial(*this) -= x; }
  Polynomial operator*(const T &x) const { return Polynomial(*this) *= x; }
  Polynomial operator/(const T &x) const { return Polynomial(*this) /= x; }
  T eval(T x) {
    T ret = 0;
    T powX = 1;
    for (long long i = 0; i < vector<T>::size(); i++) {
      ret += vector<T>::operator[](i) * powX;
      powX *= x;
    }
    return ret;
  }
};
template <class T>
Polynomial<T> simpleInterpolate(vector<T> x, vector<T> y) {
  Polynomial<T> ret;
  for (long long i = 0; i < x.size(); i++) {
    Polynomial<T> prod = {1};
    T div = 1;
    for (long long j = 0; j < x.size(); j++)
      if (i != j) {
        div *= x[i] - x[j];
        Polynomial<T> tmp = {-x[j], 1};
        prod *= tmp;
      }
    ret += prod * y[i] / div;
  }
  return ret;
}
using poly = Polynomial<mint>;
long long N, T;
long long X[55], Y[55];
mint calcNumCovered(long long t) {
  vector<long long> ys;
  for (long long i = 0; i < (N); i++) {
    ys.push_back(Y[i] - t);
    ys.push_back(Y[i] + t + 1);
  }
  sort((ys).begin(), (ys).end());
  ys.erase(unique((ys).begin(), (ys).end()), ys.end());
  vector<pair<long long, long long> > qs;
  for (long long i = 0; i < (N); i++) {
    qs.push_back({X[i] - t, i * 2});
    qs.push_back({X[i] + t + 1, i * 2 + 1});
  }
  sort((qs).begin(), (qs).end());
  long long lastX = 0;
  mint ret = 0;
  vector<long long> num(ys.size());
  for (long long i = 0; i < (qs.size()); i++) {
    long long x = qs[i].first;
    long long k = qs[i].second / 2;
    long long sum = 0;
    for (long long j = 0; j + 1 < ys.size(); j++)
      if (num[j]) sum += ys[j + 1] - ys[j];
    ret += mint(sum) * (x - lastX);
    lastX = x;
    long long l = lower_bound((ys).begin(), (ys).end(), Y[k] - t) - ys.begin();
    long long r =
        lower_bound((ys).begin(), (ys).end(), Y[k] + t + 1) - ys.begin();
    long long coef;
    if (qs[i].second & 1)
      coef = -1;
    else
      coef = 1;
    for (long long j = l; j < r; j++) num[j] += coef;
  }
  return ret;
}
vector<mint> lag(vector<mint> &x, vector<mint> &y) {
  vector<mint> ret(3);
  mint b = (x[0] - x[1]) * (x[0] - x[2]);
  b = b.inv();
  ret[0] = x[1] * x[2];
  ret[1] = -(x[1] + x[2]);
  ret[2] = 1;
  for (long long i = 0; i < (3); i++) ret[i] *= y[0] * b;
  return ret;
}
vector<mint> Lag(vector<mint> &x, vector<mint> &y) {
  vector<mint> p(3);
  for (long long j = 0; j < 3; j++) {
    vector<mint> pp = lag(x, y);
    for (long long k = 0; k < (3); k++) p[k] += pp[k];
    rotate(x.begin(), x.begin() + 1, x.end());
    rotate(y.begin(), y.begin() + 1, y.end());
  }
  return p;
}
mint calcSum(vector<mint> &p, long long l, long long r) {
  mint ret;
  ret += p[0] * (r - l);
  ret += p[1] * (mint(r - 1) * r - mint(l - 1) * l) / 2;
  ret += p[2] *
         (mint(r - 1) * r * (2 * r - 1) - mint(l - 1) * l * (2 * l - 1)) / 6;
  return ret;
}
signed main() {
  cin >> N >> T;
  for (long long i = 0; i < (N); i++) cin >> X[i] >> Y[i];
  mint ans;
  vector<long long> ts;
  for (long long i = 0; i < N; i++) {
    for (long long j = i + 1; j < N; j++) {
      long long ma = max(abs(X[i] - X[j]), abs(Y[i] - Y[j]));
      ma = (ma + 1) / 2;
      if (ma > T) continue;
      ts.push_back(ma);
    }
  }
  ts.push_back(0);
  ts.push_back(T);
  sort((ts).begin(), (ts).end());
  ts.erase(unique((ts).begin(), (ts).end()), ts.end());
  for (long long i = 0; i + 1 < ts.size(); i++) {
    long long l = ts[i];
    long long r = ts[i + 1];
    vector<mint> x(4), y(4);
    x[0] = l - 1;
    y[0] = 0;
    for (long long i = 1; i < 4; i++) {
      x[i] = l + i - 1;
      y[i] = y[i - 1] + calcNumCovered(l + i - 1);
    }
    poly p = simpleInterpolate(x, y);
    ans -= p.eval(r - 1);
  }
  ans += calcNumCovered(T) * T;
  cout << ans << endl;
}
