/**
 * DATA STRUCTURES
 */

typedef struct MM {
    float   x;
    float   y;
    int     carbone;
    int     dormance;
} MM;

typedef struct OM {
    float   x;
    float   y;
    int     carbone;
    int     lock;
} OM;

typedef struct MiorWorld {
    int     nbMM;
    int     nbOM;
    int     RA;
    float   RR;
    float   GR;
    float   K;
    int     width;
    int     minSize;
    int     CO2;
    int     lock;
} MiorWorld;


typedef struct RandomGen {
	ulong a;
	ulong b;
	ulong c;
} RandomGen;

/**
 * HELPER FUNCTIONS
 */

// Compatibility macros for OpenCL <= 1.1

#if __OPENCL_VERSION__ <= CL_VERSION_1_0
    #define atomic_inc(p) atom_inc(p)
    #define atomic_add(p, v) atom_add(p, v)
    //#define atomic_cmpxchg(p, o, v) atom_cmpxchg(p, o, v)
#endif

// End of compatibility macros

//#define LOCK_P(lock, id) while (atomic_cmpxchg(&(lock), -1, id) != id) {}
//#define LOCK_V(lock, id) while (atomic_cmpxchg(&(lock), id, -1) != -1) {}
#define store_CO2(world, val) atomic_add(&(world->CO2), val)
#define val(arr, iMM, iOM) arr[iMM * world->nbOM + iOM]

// RNG

void rng_init(RandomGen *r, long seed)
{
    r->a = seed;
    r->b = 0;
    r->c = 362436;
}

unsigned long rng_rand(RandomGen *r)
{
    const long old = r->b;
    r->b = r->a * 1103515245 + 12345;
    r->a = (~old ^ (r->b >> 3)) - r->c++;
    return r->b;
}

float rng_rand_01(RandomGen *r)
{
    return (rng_rand(r) & 4294967295) / 4294967295.0f;
}

void shuffle_indexes(global int * array, int size, RandomGen *r)
{   
    for (int i = size - 1; i >= 1; i--) {
		const int j = rng_rand_01(r) * (i + 1);
		
		const int t = array[j];
		array[j] = array[i];
		array[i] = t;
	}
}

/**
 * KERNELS
 */

kernel void topology(
    global MM *mmList,
    global OM *omList,
    global int *topo,
    global MiorWorld *world,
    global int *totals)
{
    const int iMM = get_global_id(0);
    const int iOM = get_global_id(1);
    
	const float dx = omList[iOM].x - mmList[iMM].x;
	const float dy = omList[iOM].y - mmList[iMM].y;
	const float dist = hypot(dx, dy);
	
	if (dist <= world->RA) {
	    val(topo, iMM, iOM) = 0;
	    atomic_inc(totals + iOM);
	} else {
	    val(topo, iMM, iOM) = -1;
	}
}

kernel void carbon_scatter(
    global MM * mmList,
    global OM * omList,
    global int * topo,
    global MiorWorld * world,
    global int * totals)
{
    const int iOM = get_global_id(0);
    const int total = totals[iOM];
    
    if (total == 0) {
        return;
    }
    
    const int part = (world->K * omList[iOM].carbone) / total; 
    int totalCarbon = 0;
    
    for (int iMM = 0; iMM < world->nbMM; iMM++) {
        const int offer = val(topo, iMM, iOM);
        
        if (offer >= 0) {
            val(topo, iMM, iOM) = part;
            totalCarbon += part;
        }
    }
    
    omList[iOM].carbone -= totalCarbon;
}

kernel void carbon_reduce(
    global MM * mmList,
    global OM * omList,
    global int * topo,
    global MiorWorld * world)
{
    const int iOM = get_global_id(0);
    int totalCarbon = 0;
    
    for (int iMM = 0; iMM < world->nbMM; iMM++) {
        const int part = val(topo, iMM, iOM);
        
        if (part >= 0) {
            val(topo, iMM, iOM) = 0;
            totalCarbon += part;
        }
    }
    
    omList[iOM].carbone += totalCarbon;
}

kernel void live(
    global MM *mmList,
    global OM *omList,
    global int *topo,
    global MiorWorld *world,
    global int *omIndexes,
    long   seed)
{   

    int iMM = get_global_id(0);
    global int * targets = topo + iMM * world->nbOM;
    global int * indexes = omIndexes + iMM * world->nbOM;
    RandomGen rng;
    
    for (int i = 0; i < world->nbOM; i++) {
        indexes[i] = i;
    }
    
    if (seed != 0) {
        rng_init(&rng, seed + iMM);
        shuffle_indexes(indexes, world->nbOM, &rng);
    }
    
    
    const int breathNeed = world->RR * mmList[iMM].carbone;
    int remainingBreathNeed = breathNeed;
    int iOM;
    
    // Check breath requirement
    for (int i = 0; i < world->nbOM; i++) {
        iOM = indexes[i];
        
        if (val(topo, iMM, iOM) >= 0) {
            remainingBreathNeed -= val(topo, iMM, iOM);
            
            if (remainingBreathNeed <= 0) {
                break;
            }
        }
    }
    
    // Go into dormancy if breath requirement is not met
    if (remainingBreathNeed > 0) {
        mmList[iMM].dormance = 1;
        return;
    } else {
        mmList[iMM].dormance = 0;
    }
    
    // Actual breathing and growth
    remainingBreathNeed = breathNeed;
    const int growthNeed  = world->GR * mmList[iMM].carbone + breathNeed;
    int remainingGrowthNeed = growthNeed;
    int consum = 0, growth = 0;
    
    for (int i = 0; i < world->nbOM; i++) {
        iOM = indexes[i];
        
        int offer = val(topo, iMM, iOM);
        
        if (offer >= 0) {
            // If breathing phase is not finished
            if (remainingBreathNeed > 0) {
                consum = min(offer, remainingBreathNeed);
                atomic_add(&(world->CO2), consum);
                remainingBreathNeed -= consum;
                offer -= consum;
            }
            
            consum = min(offer, remainingGrowthNeed);
            remainingGrowthNeed -= consum;
            growth += consum;
            offer -= consum;
            
            val(topo, iMM, iOM) = offer;
        }
        
        if (remainingGrowthNeed <= 0) {
            break;
        }
    }
    
    mmList[iMM].carbone += growth;
}