class Relax {
  final static int N = 6;
  final static int T = 3;
  final static double THRESHOLD = 0.1;
  final static RectDomain<2> all = [0:N-1,0:N-1];
  final static RectDomain<2> padded = [-1:N,-1:N];
  final static RectDomain<2> one = [0:T-1,0:T-1];
  final static RectDomain<2> one_padded = [-1:T,-1:T];
  final static Domain<2> xoverlap = ([-1:-1,0:T-1] + [T:T,0:T-1]
                                     + [0:T-1,-1:-1]+ [0:T-1,T:T]);
  final static Domain<2> shared = ([0:0,0:T-1] + [T-1:T-1,0:T-1]
                                   + [0:T-1,0:0]+ [0:T-1,T-1:T-1]);
  final static Domain<2> neighbors = ([-1:-1,0:0] + [0:0,1:1]
                                      + [0:0,-1:-1] + [1:1,0:0]);
  
  static void show(double [2d] m) {
    for (int j = 0; j < N; j++) {
      for (int i = 0; i < N; i++) {
        System.out.print(m[j,i] + " ");
      }
      System.out.println("");
    }
    System.out.println("");
  }
  
  public static void main(String[] args) {
    double [2d] m, my_m, mine, next;
    double single overall_max_delta;
    double max_delta;

    if (Ti.thisProc() == 0) {
      m = new double[padded];
      foreach (p in padded) {
        m[p] = 0.0;
      }
      foreach (p in [0:N-1,0:0]) {
        m[p] = 1.0;
      }
    }
    m = broadcast m from 0;

    int start_j = (Ti.thisProc() / (N/T)) * T;
    int start_i = (Ti.thisProc() % (N/T)) * T;
    my_m = m.translate([-start_j, -start_i]);

    mine = new double[one_padded];
    next = new double[one];

    mine.copy(my_m, one);

    do {
      //Ti.barrier(); // redundant
      mine.copy(my_m, xoverlap);
      Ti.barrier();

      // For testing:
      if (false) {
        my_m.copy(mine, one);
        Ti.barrier();
        if (Ti.thisProc() == 0) {
          show(m);
        }
        Ti.barrier();
      }

      max_delta = 0;

      foreach (p in one) {
        double sum = 0;
        foreach (p2 in neighbors) {
          sum += mine[p + p2];
        }
        next[p] = sum / 4;
        max_delta = Math.max(max_delta, Math.abs(next[p] - mine[p]));
      }

      mine.copy(next, one);

      my_m.copy(mine, shared);
      
      overall_max_delta = Reduce.max(max_delta);
    } while (overall_max_delta > THRESHOLD);

    my_m.copy(mine, one);
    Ti.barrier();
    if (Ti.thisProc() == 0)
      show(m);
  }
}