// Licensed under the MIT License <http://opensource.org/licenses/MIT>.
// SPDX-License-Identifier: MIT
// Copyright (c) 2021 Noah H. and Tom H.

#include <iostream>
#include <numeric>
#include <map>
#include <string>
#include <iterator>
#include <queue>
#include <random>
#include <algorithm>
#include <vector>
#include <unordered_map>
#include <functional>
#include <tuple>
#include <fstream>
#include <sched.h>
#include <stdlib.h>
#include <stdio.h>
#include <list>
#include <cstdlib>
#include <cstdio>
#include <memory>
#include <stack>
#include "robin_hood.h"
#define MS 662891



namespace hash_tuple
{
template <typename TT>
struct hash
{
    size_t operator()(TT const& tt) const
    {
        return robin_hood::hash<TT>()(tt);
    }
};
}
namespace hash_tuple{

namespace
    {
    template <class T>
    inline void hash_combine(std::size_t& seed, T const& v)
    {
        seed ^= hash_tuple::hash<T>()(v) + 0x9e3779b9 + (seed<<6) + (seed>>2);
    }
}

}

namespace hash_tuple{

    namespace
    {
        // Recursive template code derived from Matthieu M.
        template <class Tuple, size_t Index = std::tuple_size<Tuple>::value - 1>
        struct HashValueImpl
        {
          static void apply(size_t& seed, Tuple const& tuple)
          {
            HashValueImpl<Tuple, Index-1>::apply(seed, tuple);
            hash_combine(seed, std::get<Index>(tuple));
          }
        };

        template <class Tuple>
        struct HashValueImpl<Tuple,0>
        {
          static void apply(size_t& seed, Tuple const& tuple)
          {
            hash_combine(seed, std::get<0>(tuple));
          }
        };
    }

    template <typename ... TT>
    struct hash<std::tuple<TT...>>
    {
        size_t
        operator()(std::tuple<TT...> const& tt) const
        {
            size_t seed = 0;
            HashValueImpl<std::tuple<TT...> >::apply(seed, tt);
            return seed;
        }
    };

}




using namespace std;



typedef tuple<long long,long long,long long> Vertex;

void write_to_file(vector<tuple<int,int,double>> &vec, string file_name)
{
    std::ofstream f(file_name);
    for(vector<tuple<int,int,double>>::const_iterator i = vec.begin(); i != vec.end(); ++i) {
        f << get<0>(*i) << ','<< get<1>(*i) << ','<<get<2>(*i) << '\n';
    }
}



class Node
{
    public:
        Vertex v;
        Vertex temp_coords;
        //y,x,x_i,y_i
        Node *neighbours[4] = {NULL,NULL,NULL,NULL};
        Node *from = NULL;
        unsigned long long int counter = 0;
        int index=0;
        int is_leaf()
        {
            int count = -1;
            for (int ii=0;ii<4;ii++)
            {
                if (neighbours[ii]!=NULL)
                {
                    if(count!=-1)
                    {
                        return -1;
                    }
                    count=ii;
                }
                
            }

            return count;
        
        }
        vector<Vertex> coord_neighbours()
        {
            long long x = get<0>(v);
            long long y = get<1>(v);
            long long z = get<2>(v);
            vector <Vertex> ret{make_tuple(x,y+1,z+x),
                        make_tuple(x+1,y,z),
                         make_tuple(x-1,y,z),
                          make_tuple(x,y-1,z-x)};
            return ret;
        }
        Vertex coord_neighbour(int i)
        {
            long long x = get<0>(v);
            long long y = get<1>(v);
            long long z = get<2>(v);

            switch(i)
            {
                case 0: return make_tuple(x,y+1,z+x);break;
                case 1: return make_tuple(x+1,y,z);break;
                case 2: return make_tuple(x-1,y,z);break;
                default: return make_tuple(x,y-1,z-x);
            }
        }
        
};
class circular_buffer {
private:
    Node* buf_[MS];
    size_t head_ = 0;
    size_t tail_ = 0;
    const size_t max_size_ = MS;
    bool full_ = 0;
    
public:
    void put(Node* item);
    Node* get();
    void reset();

    bool empty() const;

    bool full() const;

    size_t capacity();

    size_t size() const;

};


void circular_buffer::put(Node* item)
{
    buf_[head_] = item;

    if(full_)
    {
        if(++tail_>=max_size_)
        {
            tail_ -= max_size_;
        }
    }

    if(++head_>=max_size_)
    {
        head_ -= max_size_;
    }

    full_ = head_ == tail_;
}
Node* circular_buffer::get()
{
    if(empty())
    {
        return NULL;
    }

    //Read data and advance the tail (we now have a free space)
    Node* val = buf_[tail_];
    full_ = false;
    
    if(++tail_>=max_size_)
    {
        tail_ -= max_size_;
    }

    return val;
}
void circular_buffer::reset()
{
    head_ = tail_;
    full_ = false;
}
bool circular_buffer::empty() const
{
    //if head and tail are equal, we are empty
    return (!full_ && (head_ == tail_));
}
bool circular_buffer::full() const
{
    //If tail is ahead the head by 1, we are full
    return full_;
}
size_t circular_buffer::capacity()
{
    return max_size_;
}
size_t circular_buffer::size() const
{
    size_t size = max_size_;

    if(!full_)
    {
        if(head_ >= tail_)
        {
            size = head_ - tail_;
        }
        else
        {
            size = max_size_ + head_ - tail_;
        }
    }

    return size;
}



class circular_buffer_short {
private:
    int buf_[MS];
    size_t head_ = 0;
    size_t tail_ = 0;
    const size_t max_size_ = MS;
    bool full_ = 0;
    
public:
    void put(int item);
    int get();
    void reset();

    bool empty() const;

    bool full() const;

    size_t capacity();

    size_t size() const;
    
    int operator [](size_t);

};

int circular_buffer_short::operator [](size_t i)
 {
     size_t get = i +tail_;
     if(get>=max_size_)
     {
         get -= max_size_;
     }
     return buf_[get];
 }

void circular_buffer_short::put(int item)
{
    buf_[head_] = item;

    if(full_)
    {
        if(++tail_>=max_size_)
        {
            tail_ -= max_size_;
        }
    }

    if(++head_>=max_size_)
    {
        head_ -= max_size_;
    }

    full_ = head_ == tail_;
}
int circular_buffer_short::get()
{
    if(empty())
    {
        return -1;
    }

    //Read data and advance the tail (we now have a free space)
    int val = buf_[tail_];
    full_ = false;
    
    if(++tail_>=max_size_)
    {
        tail_ -= max_size_;
    }

    return val;
}
void circular_buffer_short::reset()
{
    head_ = tail_;
    full_ = false;
}
bool circular_buffer_short::empty() const
{
    //if head and tail are equal, we are empty
    return (!full_ && (head_ == tail_));
}
bool circular_buffer_short::full() const
{
    //If tail is ahead the head by 1, we are full
    return full_;
}
size_t circular_buffer_short::capacity()
{
    return max_size_;
}
size_t circular_buffer_short::size() const
{
    size_t size = max_size_;

    if(!full_)
    {
        if(head_ >= tail_)
        {
            size = head_ - tail_;
        }
        else
        {
            size = max_size_ + head_ - tail_;
        }
    }

    return size;
}

std::function<Vertex (Vertex&)> from_to(Vertex &source, Vertex &dest,bool reflect=false)
{
    long long a1 = get<0>(source);
    long long b1 = get<1>(source);
    long long c1 = get<2>(source);

    long long a2 = get<0>(dest);
    long long b2 = get<1>(dest);
    long long c2 = get<2>(dest);
    if(reflect)
    {
        return [a1,a2,b1,b2,c1,c2](Vertex& v)
        {
            long long a3 = get<0>(v);
            long long b3 = get<1>(v);
            long long c3 = get<2>(v);
            long long aa= a2-a3;
            long long bb = b1-b3;
            return make_tuple(a1 + aa,
                              b2 +bb,
                              aa*bb + c1 + c2 - c3);
        };
    }
    else
    {
        return [a1,a2,b1,b2,c1,c2](Vertex& v)
        {
            long long a3 = get<0>(v);
            long long b3 = get<1>(v);
            long long c3 = get<2>(v);
            long long bb = b1-b3;
            long long aa = a1-a2;
            return make_tuple( a3-aa,
                               b2-bb,
                               aa*bb - c1 + c2 + c3);
        };
    }
}






class Graph
{
    public:
    vector<Node> nodes;
    unsigned long long int counter = 0;
    std::random_device rd;
    std::mt19937 rng;
    std::mt19937_64 rng_64;
    
    circular_buffer *Q1 = new circular_buffer;
    circular_buffer *Q2 = new circular_buffer;
    circular_buffer_short *Q3 = new circular_buffer_short;
    circular_buffer_short *Q4 = new circular_buffer_short;
    
    robin_hood::unordered_map <Vertex,int,hash_tuple::hash<Vertex>> coords_to_node;
    robin_hood::unordered_map <Vertex,int,hash_tuple::hash<Vertex>> kn;
    
    Graph(int size): rng(rd()),rng_64(rd())
    {
        nodes.resize(size);
        coords_to_node.reserve(size);
        
        int N = size%2==0 ? size/2 : (size+1)/2;
        int N_ = size%2==0 ? N : N-1;
        
        for (int i=0;i<N;i++)
        {
            Vertex coord = make_tuple(i,0,0);
            nodes[i].v = coord;
            nodes[i].index=i;
            coords_to_node[coord] = i;
            nodes[i].neighbours[1] = &(nodes[i+1]);
            nodes[i+1].neighbours[2] = &(nodes[i]);
        }
        Vertex coord = make_tuple(N,0,0);
        nodes[N].v = coord;
        nodes[N].index = N;
        coords_to_node[coord]=N;
        
        for (int i=1;i<N_;i++)
        {
            Vertex coord = make_tuple(i,1,i);
            nodes[N+i].v = coord;
            nodes[N+i].index = N+i;
            coords_to_node[coord] = N+i;
            nodes[N+i].neighbours[3] = &(nodes[i]);
            nodes[i].neighbours[0] = &(nodes[N+i]);
        }
        
    }

    void increment_counter(int num)
    {
        if (counter>= numeric_limits<unsigned long int>::max() - num - 1)
        {
            counter = 0;
            for(std::vector<Node>::iterator node = nodes.begin(); node != nodes.end(); ++node)
            {
                    node->counter = 0;
            }
        }
        counter = counter + num;
    }

    void run_algorithm(int num_samples,int steps_between_samples,int initial_delay,int samples_per_file,string file_name)
    {
        long long so_far = -initial_delay;
        long long num_steps=(long long)num_samples*(long long)steps_between_samples+1;
        int so_far_samples = 0;
        
        int siz = nodes.size();
        long long mult = 4*siz;
        
        std::uniform_int_distribution<long long> uni_nodes(0,mult*mult-1);

        vector<tuple<int,int,double>>ret;
        ret.resize(samples_per_file);
        int swap_count = 0;
        int cut_and_paste_count = 0;
        int s_and_i_count = 0;
        
        int swap_suc_count = 0;
        int cut_and_paste_suc_count = 0;
        
        int k=0;
        int z = 0;
        while(so_far<num_steps)
        {
            int div =so_far/(long long)steps_between_samples;
            if(so_far>=0 && div*steps_between_samples==so_far)
            {
                if(so_far>0 && k==samples_per_file)
                {
                    cout<<div<<"/"<<num_samples<<"\n";
                    write_to_file(ret, file_name+"_"+to_string(z));
                    k=0;
                    z++;
                }
                
                pair<int,int> p = longest_path();
                double b = average_branching();
                ret[k] = make_tuple(p.first,p.second,b);
                so_far_samples+=1;
                k++;
            }
            so_far+=1;

            long long random_tot = uni_nodes(rng_64);
            int random_s =random_tot/mult;
            
            int node_index = random_s%siz;
            int direction = random_s/siz;

            Node &node = nodes[node_index];
            Node *neighbour = node.neighbours[direction];
            if(neighbour != NULL)
            {
                cut_and_paste_count+=1;
                cut_and_paste_suc_count+= cut_and_paste(&node,direction);
                continue;
            }

            Vertex n_coords = node.coord_neighbour(direction);
            auto n_index_it = coords_to_node.find(n_coords);

            if(n_index_it != coords_to_node.end())
            {
                int n_index = n_index_it->second;
                Node *neighbour = &nodes[n_index];
                s_and_i_count+=1;
                search_and_insert(&node,neighbour,direction);
                continue;
            }
            swap_count+=1;
            
            
            random_s = random_tot%mult;
            int node_index_2 = random_s%siz;
            int direction_2 = random_s/siz;

            Node &node_2 = nodes[node_index_2];
            Node *node_2_neighbour = node_2.neighbours[direction_2];

            if(node_2_neighbour==NULL)
            {
                continue;
            }
            

            
            int ok_dir = node_2.is_leaf();
            if(ok_dir!=-1)
            {
                if(node_index == node_index_2)
                {
                    continue;
                }
                
                swap(node, direction, n_coords ,node_2,direction_2);
                swap_suc_count+=1;
                continue;
            }
            
            
            int ok_dir_2 = node_2_neighbour->is_leaf();
            if(ok_dir_2 == -1)
            {
                continue;
            }
            if(node_2_neighbour == &node)
            {
                continue;
            }
            swap_suc_count+=1;
            swap(node, direction, n_coords, *node_2_neighbour, ok_dir_2);
        }

        
        
    }
    
    void swap(Node &node1, int direction1,Vertex &n_coords, Node &node2, int direction2)
    {
        //cout<<"swap\n";
        Node *n2_n = node2.neighbours[direction2];
        n2_n->neighbours[3-direction2] = NULL;
        
        node1.neighbours[direction1] = &node2;
        node2.neighbours[direction2] = NULL;
        node2.neighbours[3-direction1] = &node1;
        
        coords_to_node.erase(node2.v);
        coords_to_node[n_coords]=node2.index;
        node2.v = n_coords;
        
    }

    
    bool cut_and_paste(Node *X1,int ii)
    {
        //cout<<"cut and paste\n";
        
        Node *Y1 = X1->neighbours[ii];
        
        X1->neighbours[ii] = NULL;
        Y1->neighbours[3-ii] = NULL;
        
        (*Q1).reset();
        (*Q2).reset();
        
        (*Q1).put(X1);
        (*Q2).put(Y1);
        

        int dist = 0;
        
        increment_counter(2);

        unsigned long int c1 = counter;
        unsigned long int c2 = counter-1;

        X1 -> counter = c1;
        Y1 -> counter = c2;
        
        bool XWon = false;
        
        (*Q3).reset();
        (*Q4).reset();
        while(true)
        {

            dist++;
            
            Node *X = (*Q1).get();
            Node *Y = (*Q2).get();
            

            
            if(X==NULL)
            {
                XWon=true;
                break;
            }
            if(Y==NULL)
            {
                break;
            }
            (*Q3).put(X->index);
            (*Q4).put(Y->index);
            
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = X->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c1)
                    {
                        continue;
                    }
                    neighbour->counter = c1;
                    (*Q1).put(neighbour);
                }
                
            }
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = Y->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c2)
                    {
                        continue;
                    }
                    (*Q2).put(neighbour);
                    neighbour->counter = c2;
                }
                
            }
        }
        Node* small_one;
        circular_buffer_short* small_list = NULL;
        unsigned long cc;
        if(XWon)
        {
            small_list= &(*Q3);
            cc=c1;
        }
        else
        {
            small_list= &(*Q4);
            cc=c2;
        }
        int sss = small_list->size();
        
        std::uniform_int_distribution<int> uni(0,sss*4-1);
        int r = uni(rng);
        
        int r_a = r/4;
        int dir = r-4*r_a;
        small_one = &nodes[(*small_list)[r_a]];
        
        std::uniform_int_distribution<int> uni2(0,nodes.size()-1);
        while(true)
        {
            r = uni2(rng);
            if(nodes[r].counter!=cc)
            {
                break;
            }
        }
        Node *large_one =&nodes[r];
        
        Vertex starting_coord = large_one->coord_neighbour(dir);
        
        
        if(coords_to_node.find(starting_coord)!=coords_to_node.end())
        {
            X1->neighbours[ii] = Y1;
            Y1->neighbours[3-ii] = X1;
            return false;
        }
        
        increment_counter(1);
        
        unsigned long int c = counter;
        
        (*Q1).reset();
        (*Q1).put(small_one);
        
        small_one->counter = c;
        small_one->temp_coords = starting_coord;
        
        std::function<Vertex (Vertex&)> transform = from_to(small_one->v,starting_coord,false);

        while(true)
        {
            Node *X = (*Q1).get();
            if(X==NULL)
                break;
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = X->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c)
                    {
                        continue;
                    }
                    neighbour->temp_coords = transform(neighbour->v);
                    auto found = coords_to_node.find(neighbour->temp_coords);
                    if(found!=coords_to_node.end() && nodes[found->second].counter!=cc)
                    {
                        X1->neighbours[ii] = Y1;
                        Y1->neighbours[3-ii] = X1;
                        return false;
                    }
                    
                    neighbour->counter = c;
                    (*Q1).put(neighbour);
                    
                }
            }
            
        }

        for(int i=0;i<sss;i++)
        {
            Node &node = nodes[(*small_list)[i]];
            coords_to_node.erase(node.v);
            node.v = node.temp_coords;
        }
        for(int i=0;i<sss;i++)
        {
            Node &node = nodes[(*small_list)[i]];
            coords_to_node[node.v]=node.index;
        }
        
        large_one->neighbours[dir]=small_one;
        small_one->neighbours[3-dir]=large_one;
        return true;
    }
    pair<int,int> BFS(int start_index)
    {
        Node& start = nodes[start_index];
        (*Q1).reset();
        (*Q1).put(&start);
        vector<int>  distances;
        distances.resize(nodes.size());
        
        distances[start_index] = 0;
        
        increment_counter(1);
        unsigned long long c = counter;
        start.counter = c;
        
        while(true)
        {
            Node *X = (*Q1).get();
            if(X==NULL)
            {
                break;
            }
            int X_distance = distances[X->index];
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = X->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c)
                    {
                        continue;
                    }
                    neighbour->counter = c;
                    distances[neighbour->index] =X_distance+1;
                    (*Q1).put(neighbour);
                }
            }
        }
        int ss = distances.size();
        int max_dist = 0;
        int index = start_index;
        for (int i=0;i<ss;i++)
        {
            if (distances[i]>max_dist)
            {
                max_dist = distances[i];
                index = i;
            }
        }
        return make_pair(max_dist, index);
    }
    
    pair<int,int> longest_path()
    {
        pair<int,int> ret = BFS(0);
        int start_index = ret.second;
        pair<int,int> ret_2 = BFS(start_index);
        int finish_index = ret_2.second;
        
        int distance = H_dist(nodes[start_index].v,nodes[finish_index].v);
        return make_pair(ret_2.first,distance);
    }
    double average_branching()
    {
        Node& start = nodes[0];
        
        (*Q1).reset();
        (*Q1).put(&start);
        
        deque<tuple<int,int>> vd_queue;
        stack<tuple<int,int>> vd_stack;
        
        int sss = nodes.size();
        vector<int> edge_to_branch_size(4*sss);
        
        auto pos_1 = [&sss](tuple<int,int> &pair) {return sss*get<1>(pair)+get<0>(pair) ;};
        auto pos_2 = [&sss](int a, int b) {return sss*b+a ;};
        
        increment_counter(1);
        unsigned long long c =counter;
        start.counter = c;
        while(true)
        {
            Node *X = (*Q1).get();
            if(X==NULL)
            {
                break;
            }
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = X->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c)
                    {
                        continue;
                    }
                    neighbour->counter = c;
                    vd_queue.push_front(make_tuple(X->index,i));
                    vd_queue.push_back(make_tuple(neighbour->index,3-i));
                    (*Q1).put(neighbour);
                }
            }
        }

        for (deque<tuple<int,int>>::reverse_iterator it = vd_queue.rbegin(); it != vd_queue.rend(); ++it )
        {
            vd_stack.push(*it);
            while(true)
            {
                if (vd_stack.empty())
                {
                    break;
                }
                tuple<int,int> top_vd = vd_stack.top();
                int index = pos_1(top_vd);
                if(edge_to_branch_size[index] > 0)
                {
                    vd_stack.pop();
                    continue;
                }
                
                
                Node* neighbour = nodes[get<0>(top_vd)].neighbours[get<1>(top_vd)];
                int n_index = neighbour->index;
                bool can_calculate = true;
                int sum = 1;
                for (int i=0;i<4;i++)
                {
                    if(i!=3-get<1>(top_vd) && neighbour->neighbours[i]!=NULL)
                    {
                        tuple<int,int> to_h = make_pair(n_index,i);
                        int present = edge_to_branch_size[pos_1(to_h)];
                        if(present>0)
                        {
                            sum+=present;
                        }
                        else
                        {
                            can_calculate = false;
                            vd_stack.push(std::move(to_h));
                        }
                    }
                }
                if(can_calculate)
                {
                    edge_to_branch_size[index] = sum;
                    vd_stack.pop();
                }

            }
        }
        int sum = 0;
        for (int i=0;i<nodes.size()-1;i++)
        {
            tuple<int,int> top = std::move(vd_queue.front());
            int n_index =nodes[get<0>(top)].neighbours[get<1>(top)]->index;
            sum+= -1+min(edge_to_branch_size[pos_1(top)],edge_to_branch_size[pos_2(n_index,3-get<1>(top))]);
            vd_queue.pop_front();
        }
        double average = (double)sum/(double)(nodes.size()-1);
        return average;
    }
    
    
    
    int o_dist(Vertex &v)
    {
        long long x = get<0>(v);
        long long y = get<1>(v);
        long long z = get<2>(v);
        
        if (x==0 && y==0 && z==0)
        {
            return 0;
        }
        if(z<=0)
        {
            y=-y;
            z=-z;
        }
        if (abs(y)>abs(x))
        {
            long long temp = y;
            y=x;
            x=temp;
        }
        if (x<=0)
        {
            x=-x;
            y=-y;
        }

        if(y>=0)
        {
            double sz  = double(sqrt(z));
            if(x<=sz)
            {
                return 2*(ceil(2*sz))-x-y;
            }
            else
            {
                if(x*y>=z)
                {
                    return x+y;
                }
                else
                {
                    return 2*ceil(double(z)/double(x))+x-y;
                }

            }

        }
        else
        {
            double sz = sqrt(double(z-x*y));
            if (x<=sz)
            {
                return 2*ceil(2*sz)-x+y;
            }
            else
            {
                return 2*ceil(double(z)/double(x))+x-y;
            }

        }

    }
    int H_dist(Vertex v1, Vertex v2)
    {
        Vertex origin = make_tuple(0,0,0);
        Vertex translated = from_to(v1, origin)(v2);
        return o_dist(translated);
    }
    
    void search_and_insert(Node *X1, Node *Y1, int ii)
    {
        //cout<<"s_ans_i\n";
        (*Q1).reset();
        (*Q2).reset();
        
        (*Q1).put(X1);
        (*Q2).put(Y1);
        increment_counter(2);

        unsigned long int c1 = counter;
        unsigned long int c2 = counter-1;

        bool success = false;
        
        int dist = 0;

        Node *success_X = NULL;
        Node *success_Y = NULL;

        int X_d =0;
        int Y_d =0;
        
        X1 -> counter = c1;
        Y1 -> counter = c2;
        
        while(success == false)
        {
            dist++;
            
            Node *X = (*Q1).get();
            Node *Y = (*Q2).get();
            
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = X->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c1)
                    {
                        continue;
                    }
                        if (neighbour->counter == c2)
                    {
                        success_X = X;
                        success_Y = neighbour;

                        X_d = i;
                        Y_d = 3-i;

                        success = true;
                        break;
                    }
                    (*Q1).put(neighbour);
                    neighbour->from = X;
                    neighbour->counter = c1;
                }
                
            }
            if(success==true)
            {
                break;
            }
            for (int i = 0; i<4; i++)
            {
                Node *neighbour  = Y->neighbours[i];
                if(neighbour!=NULL)
                {
                    if (neighbour->counter == c2)
                    {
                        continue;
                    }
                    if (neighbour->counter == c1)
                    {
                        success_Y = Y;
                        success_X = neighbour;

                        Y_d = i;
                        X_d = 3-i;

                        success = true;
                        break;
                    }
                    (*Q2).put(neighbour);
                    neighbour->from = Y;
                    neighbour->counter = c2;
                }
                
            }
        }
        vector <Node *> hist = {};
        hist.reserve(2*dist+2);
        Node *last = success_X;

        hist.push_back(X1);
        while(last != X1)
        {
            hist.push_back(last);
            last = last->from;
        }
        
        last = success_Y;
        
        while(last != Y1)
        {
            hist.push_back(last);
            last = last->from;
        }
        std::uniform_int_distribution<int> uni(0,hist.size()-1);
        int r = uni(rng);
        if (r == 0)
        {
            success_X->neighbours[X_d] = NULL;
            success_Y->neighbours[Y_d] = NULL;
        }
        else
        {
            Node *v1 = hist[r];
            Node *v2 = hist[r]->from;
            for (int i=0;i<4;i++)
            {
                if (v1->neighbours[i] == v2)
                {
                    v1->neighbours[i] = NULL;
                    v2->neighbours[3-i] = NULL;
                    break;
                }
            }
        
        }
        X1->neighbours[ii] = Y1;
        Y1->neighbours[3-ii] = X1;
    }
    

};


// RUN LATTICE TREE SIMULATION
// 1. graph_size = number of vertices in tree
// 2. steps_between_samples_to_gs = (number of steps between samples)/(graph size)
// 3. num_samples = number of samples
// 4. samples_per_file = samples per file
// 5. file_name = output path

int main(int argc, char** argv)
{
    cout<<"start\n";
    
    int graph_size = stoi(argv[1]);
    if (graph_size>MS)
    {
        fprintf(stderr, "Graph size too large.\n");
        return -1;
    }
    int steps_between_samples_to_gs = stoi(argv[2]);
    
    int steps_between_samples = graph_size*steps_between_samples_to_gs;
    
    int num_samples = stoi(argv[3]);

    int samples_per_file = stoi(argv[4]);
    
    string file_name = argv[5];
    
    Graph G(graph_size);

    
G.run_algorithm(num_samples,steps_between_samples,2*steps_between_samples,samples_per_file,file_name);
    cout<<"done\n";
}