Language: C++

#include <bits/stdc++.h>
using namespace std;
#define _MACRO(_1, _2, _3, NAME, ...) NAME
#define _repl(i,a,b) for(int i=(int)(a);i<(int)(b);i++)
#define _rep(i,n) _repl(i,0,n)
#define rep(...) _MACRO(__VA_ARGS__, _repl, _rep)(__VA_ARGS__)
#define pb push_back
#define all(x) begin(x),end(x)
#define uniq(x) sort(all(x)),(x).erase(unique(all(x)),end(x))
#ifdef LOCAL
#define dbg(...) _dbg(#__VA_ARGS__, __VA_ARGS__)
void _dbg(string){cerr<<endl;}
template<class H,class... T> void _dbg(string s,H h,T... t){int l=s.find(',');cerr<<s.substr(0,l)<<" = "<<h<<", ";_dbg(s.substr(l+1),t...);}
template<class T,class U> ostream& operator<<(ostream &o, const pair<T,U> &p){o<<"("<<p.first<<","<<p.second<<")";return o;}
template<class T> ostream& operator<<(ostream &o, const vector<T> &v){o<<"[";for(T t:v){o<<t<<",";}o<<"]";return o;}
constexpr int SIZE = 100;
#else
#define dbg(...) {}
constexpr int SIZE = 4096;
#endif

int val[2][2005][2005];

int main(){
  int n,d;
  cin>>n>>d;
  vector<vector<int>> v(d,vector<int>(d,0));
  rep(i,n){
    int x,y;
    cin>>x>>y;
    v[x%d][y%d]++;
  }

  long l = 1, r = n*d;
  while(r-l > 1) {
    long mid = (l+r)/2;
    auto f = [&](){
      if (mid * mid < n) return false;
      rep(i,d) rep(j,d) if(v[i][j] > mid/d * mid/d) return false;
      if(mid%d == 0) return true;

      fill(val[0], val[2], 0);
      rep(i,1,d+1) rep(j,1,d+1){
        val[0][i][j] = val[0][i+d][j] = val[0][i][j+d] = val[0][i+d][j+d] = v[i-1][j-1] <= (mid/d + 1)*(mid/d);
        val[1][i][j] = val[1][i+d][j] = val[1][i][j+d] = val[1][i+d][j+d] = v[i-1][j-1] <= (mid/d + 1)*(mid/d + 1);
      }
      rep(k,2){
        rep(i,2*d+1){
          rep(j,2*d+1){
            if(i>0) val[k][i][j] += val[k][i-1][j];
            if(j>0) val[k][i][j] += val[k][i][j-1];
            if(i>0 && j>0) val[k][i][j] -= val[k][i-1][j-1];
          }
        }
      }

      rep(i,d) rep(j,d){
        if(val[1][i + mid%d][j + mid%d] - val[1][i][j + mid%d] - val[1][i + mid%d][j] + val[1][i][j] < (mid%d)*(mid%d)) continue;
        if(val[0][i + mid%d][j + d] - val[0][i][j + d] - val[0][i + mid%d][j + mid%d] + val[0][i][j + mid%d] < (d - mid%d)*(mid%d)) continue;
        if(val[0][i + d][j + mid%d] - val[0][i + d][j] - val[0][i + mid%d][j + mid%d] + val[1][i + mid%d][j] < (d - mid%d)*(mid%d)) continue;
        return true;
      }
      return false;
    };
    if (f()) r = mid;
    else l = mid;
  }

  cout << l << endl;

  return 0;
}


Compiler output:
In file included from /usr/include/c++/14/algorithm:60,
                 from /usr/include/x86_64-linux-gnu/c++/14/bits/stdc++.h:51,
                 from a.cc:1:
/usr/include/c++/14/bits/stl_algobase.h: In instantiation of 'typename __gnu_cxx::__enable_if<std::__is_scalar<_Tp>::__value, void>::__type std::__fill_a1(_ForwardIterator, _ForwardIterator, const _Tp&) [with _ForwardIterator = int (*)[2005]; _Tp = int; typename __gnu_cxx::__enable_if<__is_scalar<_Tp>::__value, void>::__type = void]':
/usr/include/c++/14/bits/stl_algobase.h:998:21:   required from 'void std::__fill_a(_FIte, _FIte, const _Tp&) [with _FIte = int (*)[2005]; _Tp = int]'
  998 |     { std::__fill_a1(__first, __last, __value); }
      |       ~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~
/usr/include/c++/14/bits/stl_algobase.h:1029:20:   required from 'void std::fill(_ForwardIterator, _ForwardIterator, const _Tp&) [with _ForwardIterator = int (*)[2005]; _Tp = int]'
 1029 |       std::__fill_a(__first, __last, __value);
      |       ~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~
a.cc:42:11:   required from here
   42 |       fill(val[0], val[2], 0);
      |       ~~~~^~~~~~~~~~~~~~~~~~~
/usr/include/c++/14/bits/stl_algobase.h:952:18: error: incompatible types in assignment of 'const int' to 'int [2005]'
  952 |         *__first = __tmp;
      |         ~~~~~~~~~^~~~~~~