particles/C/terrain.h

#ifndef TERRAIN
#define TERRAIN

#include <stdlib.h>
#include <omp.h>

typedef struct world_coord_t {
    size_t i;
    size_t j;
} world_coord_t;

typedef struct terrain_t {
    size_t pop;
    float pheromone;
    float decay;
    float increase;
} terrain_t;

void terrain_update(terrain_t *terrain, size_t width, size_t height, float time_delta){
    // Extra decay on borders
    #pragma omp parallel for
    for(size_t i = 0 ; i < height ; ++i){
        terrain[i * width].pheromone -= 2 * terrain[i * width].decay * time_delta;
        terrain[(i + 1) * width - 1].pheromone -= 2 * terrain[(i + 1) * width - 1].decay * time_delta;
    }
    #pragma omp parallel for
    for(size_t j = 1 ; j < width - 1 ; ++j){
        terrain[j].pheromone -= 2 * terrain[j].decay * time_delta;
        terrain[(height - 1) * width + j].pheromone -= 2 * terrain[(height - 1) * width + j].decay * time_delta;
    }

    // Decay and increase
    #pragma omp parallel for
    for(size_t i = 0 ; i < height ; ++i){
        for(size_t j = 0 ; j < width ; ++j){
            size_t k = i * width + j;
            terrain[k].pheromone += (terrain[k].increase - terrain[k].decay) * time_delta;
            if(terrain[k].pheromone < 0.0f)
                terrain[k].pheromone = 0.0f;
        }
    }
    // diffusion
    const float pourcentage = 2.5f;
    for(size_t i = 1 ; i < height -1 ; ++i){
        for(size_t j = 1 ; j < width -1 ; ++j){
            if(terrain[i * width + j].pheromone == 0.0f)
                continue;
            
            terrain[(i - 1) * width + j].pheromone += terrain[i * width + j].pheromone * pourcentage * time_delta;
            terrain[(i + 1) * width + j].pheromone += terrain[i * width + j].pheromone * pourcentage * time_delta;
            terrain[i * width + j - 1].pheromone += terrain[i * width + j].pheromone * pourcentage * time_delta;
            terrain[i * width + j + 1].pheromone += terrain[i * width + j].pheromone * pourcentage * time_delta;
            terrain[i * width + j].pheromone -= terrain[i * width + j].pheromone * 4 * pourcentage * time_delta;
        }
    }
}

#endif //TERRAIN