diff --git a/MPI/Matrix.h b/MPI/Matrix.h
index 6154279..2bbe6c6 100644
--- a/MPI/Matrix.h
+++ b/MPI/Matrix.h
@@ -1,181 +1,188 @@
 //
 // Created by Joachim Jacques Koerfer on 06.05.2019.
 //
 
 #ifndef SEQUENTIAL_MATRIX_H
 #define SEQUENTIAL_MATRIX_H
 
 #include <string>
 #include <mpi.h>
 #include <fstream>
 #include <iostream>
 #include <chrono>
 
 #include "exception.h"
 
 /**
  * This matrix class is just a simple wrapper that tries to be as simple as possible without performances sacrifices.
  * It stores the element in a row-first format
  */
 class Matrix {
 private:
     using iterator_ = double*;
 
     int size_x_;
     int size_y_;
     double *data_;
 
     //For MPI
     int prank_;
     int psize_;
     int read_offset_y_;
     bool init_ = false;
 public:
     /**
-     * TODO : avoid empty declaration
-     */
-    Matrix() {};
-    /**
-     * Constructor
-     * @param size_x
-     * @param size_y
+     * Constructs an empty Matrix of size size_x*size_y
      * @param init if init=false, no value in the matrix is initialized
      */
     Matrix(int size_x, int size_y, bool init = true) : size_x_(size_x), size_y_(size_y), init_(true) {
         //Retrieve info of MPI
         MPI_Comm_rank(MPI_COMM_WORLD, &prank_);
         MPI_Comm_size(MPI_COMM_WORLD, &psize_);
 
         data_ = new double[size_x*size_y];
         if(init) {
             for(int i = 0;i < size_x*size_y;++i) {
                 data_[i] = 0;
             }
         }
     };
     /**
-     * Construct Matrix from a binary file
+     * Construct Matrix from a binary file with MPI-IO
      * @param file_name path to the file
      * @param size_x
      * @param size_y
      * @param init if init=false, no value in the matrix is initialized
      */
     Matrix(const std::string &file_name, int size_x) : size_x_(size_x), init_(true) {
         //Retrieve info of MPI
         MPI_Comm_rank(MPI_COMM_WORLD, &prank_);
         MPI_Comm_size(MPI_COMM_WORLD, &psize_);
         MPI_Status status;
 
         //Find the y size of buffer to be read (depends on prank)
         int n = size_x_ / psize_;
         int rest = size_x_ % psize_;
         int local_size_y = n + (prank_ < rest ? 1 : 0);
 
         size_y_ = local_size_y;
         //Add the ghost cells
         if(psize_ > 1)
             size_y_ += (prank_ > 0 && prank_ < psize_ - 1 ? 2 : 1);
 
         //Find the offset at which the file should be read
         read_offset_y_ = (prank_ < rest ? (n+1)*prank_ : rest + n*prank_);
 
         data_ = new double[size_x_*size_y_];
         MPI_File file;
         int file_open_error = MPI_File_open(MPI_COMM_WORLD, file_name.c_str(), MPI_MODE_RDONLY, MPI_INFO_NULL, &file);
 
         if (file_open_error != MPI_SUCCESS && prank_ == 0) {
             throw FileOpenError(file_open_error, file_name);
         }
 
 
         MPI_File_set_view(file, read_offset_y_*size_x_, MPI_DOUBLE, MPI_CHAR, "native", MPI_INFO_NULL);
         file_open_error = MPI_File_read_ordered(file, data_ + (prank_ != 0 ? size_x_ : 0), local_size_y*size_x_, MPI_DOUBLE, &status);
 
         if (file_open_error != MPI_SUCCESS && prank_ == 0) {
             throw FileOpenError(file_open_error, file_name);
         }
         MPI_File_close(&file);
     }
 
-    /**
-     * Copy constructor
+    /** Copy constructor
+     * Makes a copy from the parameter copy
      * @param copy
      */
     Matrix(const Matrix &copy) {
         prank_ = copy.prank_;
         psize_ = copy.psize_;
         read_offset_y_ = copy.read_offset_y_;
 
         size_x_ = copy.size_x_;
         size_y_ = copy.size_y_;
         init_ = true;
         data_ = new double[size_x_*size_y_];
         for(int i = 0;i <size_x_*size_y_;i++) {
             data_[i] = copy.data_[i];
         }
     }
+    /**
+     * Destructor
+     */
     ~Matrix() {
         if(init_)
             delete[] data_ ;
     };
 
+    /**
+     * Saves the matrix to the desired file with MPI-IO
+     * @param file_name
+     */
     void save(const std::string &file_name) {
         MPI_File file;
         MPI_Status status;
 
         MPI_File_open(MPI_COMM_WORLD, file_name.c_str(), MPI_MODE_WRONLY | MPI_MODE_CREATE, MPI_INFO_NULL, &file);
         MPI_File_set_size(file, size_x_*size_x_);
 
         MPI_File_set_view(file, (read_offset_y_+(prank_ != 0 ? 1 : 0))*size_x_, MPI_DOUBLE, MPI_CHAR, "native", MPI_INFO_NULL);
 
         int num_ghost_cells = (prank_ > 0 && prank_ < psize_ - 1 ? 2 : 1);
         if(psize_ == 1)
             num_ghost_cells = 0;
 
         std::cout << "rank " << prank_ << " ghost cells : " << num_ghost_cells << " size_y : " << size_y_ << std::endl;
 
         int file_write_error = MPI_File_write_ordered(file, (data_ + (prank_ == 0 ? 0 : size_x_)),
                 (size_y_ - num_ghost_cells)*size_x_, MPI_DOUBLE, &status);
 
         if (file_write_error != MPI_SUCCESS && prank_ == 0) {
             throw FileWriteError(file_write_error, file_name);
         }
         MPI_File_close(&file);
     }
 
     /**
      * Gets an element from the matrix by reference
      * This function doesn't verify if the indices are out of bound.
      * @param x index x
      * @param y index y
      * @return reference to the element
      */
     double& operator()(int x, int y) { return data_[size_x_*y + x];}
 
+    /**
+     * Gets the pointer to an element
+     * @param x
+     * @param y
+     * @return point to the element
+     */
     double* at(int x, int y)  { return &data_[size_x_*y + x]; }
 
     inline int get_size_x() { return size_x_; }
     inline int get_size_y() { return size_y_; }
 
     /**
      * Gets an element in the array of the matrix by reference
      * This function doesn't verify if the indices are out of bound.
      * @param index index which range from 0 to size_x*size_y
      * @return reference to the element
      */
     double& operator[](int index) { return data_[index]; };
 
     /**
      * Iterator function that is there for usage of for(auto d : Matrix) {}
      * @return reference to the beginning of the array of doubles
      */
     iterator_ begin() { return &data_[0]; }
     /**
      * Iterator function that is there for usage of for(auto px : Matrix) {}
      * @return reference to the end of the array of doubles
      */
     iterator_ end() { return &data_[size_x_*size_y_]; };
 };
 
 
 #endif //SEQUENTIAL_MATRIX_H
diff --git a/MPI/Simulation.h b/MPI/Simulation.h
index e07f865..8fac77c 100644
--- a/MPI/Simulation.h
+++ b/MPI/Simulation.h
@@ -1,84 +1,132 @@
 //
 // Created by Joachim Koerfer on 10.05.2019.
 //
 
 #ifndef SEQUENTIAL_SIMULATION_H
 #define SEQUENTIAL_SIMULATION_H
 
 #include <memory>
 #include <mpi.h>
 #include <string>
 #include "Matrix.h"
 
 #define SQRT_2 1.414213562373095
 
 struct State {
     int niter = 0;
     double T = 0;
     double max_nu = 0;
     double delta_t = 0;
 
 };
 
 typedef std::chrono::high_resolution_clock Clock;
 
 class Simulation {
 private:
     std::unique_ptr<Matrix> H_;
     std::unique_ptr<Matrix> Hu_;
     std::unique_ptr<Matrix> Hv_;
 
     std::unique_ptr<Matrix> Zdx_;
     std::unique_ptr<Matrix> Zdy_;
 
     //Temp storage
     std::unique_ptr<Matrix> Ht_;
     std::unique_ptr<Matrix> Hut_;
     std::unique_ptr<Matrix> Hvt_;
 
     //Variables for MPI
     int prank_;
     int psize_;
     int size_y_;
     int previous_proc_;
     int next_proc_;
     double communication_time_;
 
     int size_;
     int nx_;
     double delta_nx_;
     double Tend_;
     const double g_ = 127267.2;
     bool verbose_;
 
     //Variables relative to the state of the simulation
     double T_;
     double delta_t_;
     double cst_;
     int niter_;
 
     //Variables relative to the measure of performance
     Clock::time_point start_time_;
     Clock::time_point end_time_;
     double measured_time = 0;
     double nflops_ = 0;
 
+    /**
+     * Exchanges the ghost cells to its neighbors with MPI communications
+     * @param mat
+     */
     void exchange_ghost_cells(Matrix &mat);
+    /**
+     * Computes the Delta t in the numerical scheme
+     * @return delta_t
+     */
     double compute_delta_t();
+    /**
+     * Updates H, Hu and Hv according to the numerical scheme
+     */
     void update_matrices();
+    /**
+     * Apply the conditions that no cell in H must go below 0, and deactivate cells that are below a certain threshold in Hu and Hv
+     */
     void impose_tolerances();
+    /**
+     * Enforces the boundary conditions of the numerical scheme
+     */
     void enforce_boundary_conditions();
+    /**
+     * Copies the old border of H to Ht before H and Ht are swapped by the pointers (same for Hu and Hv)
+     */
     void copy_old_borders();
 
 
 public:
+    /**
+     * Initialize the simulation with the initial conditions present in predefined files
+     * @param path path to the folder of the initial conditions (H, Hu, Ht, Zdx and Zdy)
+     * @param nx size of one side of the grid
+     * @param size size in km of the problem
+     * @param Tend end time for the simulation
+     * @param verbose if true, logs the performances of the simulation in the console
+     * @param random_values if true, initializes the values to a subset of a realistic 8000x8000 grid
+     */
     Simulation(std::string path, int nx, int size, double Tend, bool verbose=true, bool random_values = false);
+    /**
+     * Computes the whole simulation
+     */
     void compute();
+    /**
+     * Computes only one step of the simulation. Return true as long as the simulation is not ended
+     * @return
+     */
     bool compute_step();
+    /**
+     * Gets the flops/s performance of the simulation (only valid after the simulation has ended)
+     * @return flops/s performance
+     */
     double get_flops_performance() { return nflops_; }
+    /**
+     * Gets the time-to-solution for the simulation
+     * @return time to solution
+     */
     double get_TTS() { return measured_time; }
+    /**
+     * Saves the simulation to the disk
+     * @param file_name
+     */
     void save(const std::string &file_name);
 };
 
 
 #endif //SEQUENTIAL_SIMULATION_H