diff --git a/modules/wfc/array_2d.h b/modules/wfc/array_2d.h
index 79b14d4aa..749c2644a 100644
--- a/modules/wfc/array_2d.h
+++ b/modules/wfc/array_2d.h
@@ -2,7 +2,7 @@
 #define FAST_WFC_UTILS_ARRAY2D_HPP_
 
 #include "assert.h"
-#include <vector>
+#include "core/vector.h"
 
 template <typename T>
 class Array2D {
@@ -10,7 +10,7 @@ public:
 	std::size_t height;
 	std::size_t width;
 
-	std::vector<T> data;
+	Vector<T> data;
 
 	Array2D(std::size_t p_height, std::size_t p_width) {
 		height = p_height;
@@ -21,7 +21,8 @@ public:
 	Array2D(std::size_t p_height, std::size_t p_width, T p_value) {
 		height = p_height;
 		width = p_width;
-		data.resize(width * height, p_value);
+		data.resize(width * height);
+		data.fill(p_value);
 	}
 
 	const T &get(std::size_t i, std::size_t j) const {
@@ -31,7 +32,7 @@ public:
 
 	T &get(std::size_t i, std::size_t j) {
 		assert(i < height && j < width);
-		return data[j + i * width];
+		return data.write[j + i * width];
 	}
 
 	Array2D<T> reflected() const {
@@ -78,18 +79,4 @@ public:
 	}
 };
 
-namespace std {
-template <typename T>
-class hash<Array2D<T>> {
-public:
-	std::size_t operator()(const Array2D<T> &a) const {
-		std::size_t seed = a.data.size();
-		for (const T &i : a.data) {
-			seed ^= hash<T>()(i) + (std::size_t)0x9e3779b9 + (seed << 6) + (seed >> 2);
-		}
-		return seed;
-	}
-};
-} //namespace std
-
 #endif
diff --git a/modules/wfc/array_3d.h b/modules/wfc/array_3d.h
index d91ead2ad..c3bb8ff6b 100644
--- a/modules/wfc/array_3d.h
+++ b/modules/wfc/array_3d.h
@@ -2,7 +2,7 @@
 #define FAST_WFC_UTILS_ARRAY3D_HPP_
 
 #include "assert.h"
-#include <vector>
+#include "core/vector.h"
 
 template <typename T>
 class Array3D {
@@ -11,7 +11,7 @@ public:
 	std::size_t width;
 	std::size_t depth;
 
-	std::vector<T> data;
+	Vector<T> data;
 
 	Array3D(std::size_t p_height, std::size_t p_width, std::size_t p_depth) {
 		height = p_height;
@@ -24,7 +24,8 @@ public:
 		height = p_height;
 		width = p_width;
 		depth = p_depth;
-		data.resize(width * height * depth, value);
+		data.resize(width * height * depth);
+		data.fill(value);
 	}
 
 	const T &get(std::size_t i, std::size_t j, std::size_t k) const {
@@ -33,7 +34,7 @@ public:
 	}
 
 	T &get(std::size_t i, std::size_t j, std::size_t k) {
-		return data[i * width * depth + j * depth + k];
+		return data.write[i * width * depth + j * depth + k];
 	}
 
 	bool operator==(const Array3D &a) const {
diff --git a/modules/wfc/overlapping_wfc.h b/modules/wfc/overlapping_wfc.h
index 414439e46..1c2b08a55 100644
--- a/modules/wfc/overlapping_wfc.h
+++ b/modules/wfc/overlapping_wfc.h
@@ -3,7 +3,7 @@
 
 #include <algorithm>
 #include <unordered_map>
-#include <vector>
+#include "core/vector.h"
 
 #include "array_2d.h"
 #include "wfc.h"
@@ -35,15 +35,15 @@ private:
 
 	OverlappingWFCOptions options;
 
-	std::vector<Array2D<T>> patterns;
+	Vector<Array2D<T>> patterns;
 
 	WFC wfc;
 
 	OverlappingWFC(
 			const Array2D<T> &input, const OverlappingWFCOptions &options,
 			const int &seed,
-			const std::pair<std::vector<Array2D<T>>, std::vector<double>> &patterns,
-			const std::vector<std::array<std::vector<unsigned>, 4>> &propagator) :
+			const std::pair<Vector<Array2D<T>>, Vector<double>> &patterns,
+			const Vector<std::array<Vector<unsigned>, 4>> &propagator) :
 			input(input), options(options), patterns(patterns.first), wfc(options.periodic_output, seed, patterns.second, propagator, options.get_wave_height(), options.get_wave_width()) {
 		// If necessary, the ground is set.
 		if (options.ground) {
@@ -53,12 +53,12 @@ private:
 
 	OverlappingWFC(const Array2D<T> &input, const OverlappingWFCOptions &options,
 			const int &seed,
-			const std::pair<std::vector<Array2D<T>>, std::vector<double>>
+			const std::pair<Vector<Array2D<T>>, Vector<double>>
 					&patterns) :
 			OverlappingWFC(input, options, seed, patterns,
 					generate_compatible(patterns.first)) {}
 
-	void init_ground(WFC &wfc, const Array2D<T> &input, const std::vector<Array2D<T>> &patterns, const OverlappingWFCOptions &options) {
+	void init_ground(WFC &wfc, const Array2D<T> &input, const Vector<Array2D<T>> &patterns, const OverlappingWFCOptions &options) {
 		unsigned ground_pattern_id = get_ground_pattern_id(input, patterns, options);
 
 		for (unsigned j = 0; j < options.get_wave_width(); j++) {
@@ -74,7 +74,7 @@ private:
 		wfc.propagate();
 	}
 
-	static unsigned get_ground_pattern_id(const Array2D<T> &input, const std::vector<Array2D<T>> &patterns, const OverlappingWFCOptions &options) {
+	static unsigned get_ground_pattern_id(const Array2D<T> &input, const Vector<Array2D<T>> &patterns, const OverlappingWFCOptions &options) {
 		// Get the pattern.
 		Array2D<T> ground_pattern = input.get_sub_array(input.height - 1, input.width / 2, options.pattern_size, options.pattern_size);
 
@@ -91,14 +91,14 @@ private:
 	}
 
 	//Return the list of patterns, as well as their probabilities of apparition.
-	static std::pair<std::vector<Array2D<T>>, std::vector<double>> get_patterns(const Array2D<T> &input, const OverlappingWFCOptions &options) {
+	static std::pair<Vector<Array2D<T>>, Vector<double>> get_patterns(const Array2D<T> &input, const OverlappingWFCOptions &options) {
 		std::unordered_map<Array2D<T>, unsigned> patterns_id;
-		std::vector<Array2D<T>> patterns;
+		Vector<Array2D<T>> patterns;
 
 		// The number of time a pattern is seen in the input image.
-		std::vector<double> patterns_weight;
+		Vector<double> patterns_weight;
 
-		std::vector<Array2D<T>> symmetries(
+		Vector<Array2D<T>> symmetries(
 				8, Array2D<T>(options.pattern_size, options.pattern_size));
 		unsigned max_i = options.periodic_input
 				? input.height
@@ -166,15 +166,14 @@ private:
 	// If agrees(pattern1, pattern2, dy, dx), then compatible[pattern1][direction]
 	// contains pattern2, where direction is the direction defined by (dy, dx)
 	// (see direction.hpp).
-	static std::vector<std::array<std::vector<unsigned>, 4>> generate_compatible(const std::vector<Array2D<T>> &patterns) {
-		std::vector<std::array<std::vector<unsigned>, 4>> compatible = std::vector<std::array<std::vector<unsigned>, 4>>(patterns.size());
+	static Vector<std::array<Vector<unsigned>, 4>> generate_compatible(const Vector<Array2D<T>> &patterns) {
+		Vector<std::array<Vector<unsigned>, 4>> compatible = Vector<std::array<Vector<unsigned>, 4>>(patterns.size());
 
 		// Iterate on every dy, dx, pattern1 and pattern2
 		for (unsigned pattern1 = 0; pattern1 < patterns.size(); pattern1++) {
 			for (unsigned direction = 0; direction < 4; direction++) {
 				for (unsigned pattern2 = 0; pattern2 < patterns.size(); pattern2++) {
-					if (agrees(patterns[pattern1], patterns[pattern2],
-								directions_y[direction], directions_x[direction])) {
+					if (agrees(patterns[pattern1], patterns[pattern2], directions_y[direction], directions_x[direction])) {
 						compatible[pattern1][direction].push_back(pattern2);
 					}
 				}
@@ -200,27 +199,27 @@ private:
 					output.get(y, x) = patterns[output_patterns.get(y, x)].get(0, 0);
 				}
 			}
+
 			for (unsigned y = 0; y < options.get_wave_height(); y++) {
-				const Array2D<T> &pattern =
-						patterns[output_patterns.get(y, options.get_wave_width() - 1)];
+				const Array2D<T> &pattern = patterns[output_patterns.get(y, options.get_wave_width() - 1)];
 				for (unsigned dx = 1; dx < options.pattern_size; dx++) {
 					output.get(y, options.get_wave_width() - 1 + dx) = pattern.get(0, dx);
 				}
 			}
+
 			for (unsigned x = 0; x < options.get_wave_width(); x++) {
-				const Array2D<T> &pattern =
-						patterns[output_patterns.get(options.get_wave_height() - 1, x)];
+				const Array2D<T> &pattern = patterns[output_patterns.get(options.get_wave_height() - 1, x)];
 				for (unsigned dy = 1; dy < options.pattern_size; dy++) {
 					output.get(options.get_wave_height() - 1 + dy, x) =
 							pattern.get(dy, 0);
 				}
 			}
-			const Array2D<T> &pattern = patterns[output_patterns.get(
-					options.get_wave_height() - 1, options.get_wave_width() - 1)];
+
+			const Array2D<T> &pattern = patterns[output_patterns.get(options.get_wave_height() - 1, options.get_wave_width() - 1)];
+
 			for (unsigned dy = 1; dy < options.pattern_size; dy++) {
 				for (unsigned dx = 1; dx < options.pattern_size; dx++) {
-					output.get(options.get_wave_height() - 1 + dy,
-							options.get_wave_width() - 1 + dx) = pattern.get(dy, dx);
+					output.get(options.get_wave_height() - 1 + dy, options.get_wave_width() - 1 + dx) = pattern.get(dy, dx);
 				}
 			}
 		}
diff --git a/modules/wfc/propagator.cpp b/modules/wfc/propagator.cpp
index 4e0d70ac6..04e6c02ba 100644
--- a/modules/wfc/propagator.cpp
+++ b/modules/wfc/propagator.cpp
@@ -23,8 +23,8 @@ void Propagator::propagate(Wave &wave) {
 	while (propagating.size() != 0) {
 		// The cell and pattern that has been set to false.
 		unsigned y1, x1, pattern;
-		std::tie(y1, x1, pattern) = propagating.back();
-		propagating.pop_back();
+		std::tie(y1, x1, pattern) = propagating[propagating.size() - 1];
+		propagating.resize(propagating.size() - 1);
 
 		// We propagate the information in all 4 directions.
 		for (unsigned direction = 0; direction < 4; direction++) {
@@ -48,24 +48,25 @@ void Propagator::propagate(Wave &wave) {
 
 			// The index of the second cell, and the patterns compatible
 			unsigned i2 = x2 + y2 * wave.width;
-			const std::vector<unsigned> &patterns =
-					propagator_state[pattern][direction];
+			const Vector<unsigned> &patterns = propagator_state[pattern][direction];
 
 			// For every pattern that could be placed in that cell without being in
 			// contradiction with pattern1
-			for (auto it = patterns.begin(), it_end = patterns.end(); it < it_end;
-					++it) {
+			int size = patterns.size();
+			for (int i = 0; i < size; ++i) {
+				unsigned int pattern = patterns[i];
+				
 				// We decrease the number of compatible patterns in the opposite
 				// direction If the pattern was discarded from the wave, the element
 				// is still negative, which is not a problem
-				std::array<int, 4> &value = compatible.get(y2, x2, *it);
+				std::array<int, 4> &value = compatible.get(y2, x2, pattern);
 				value[direction]--;
 
 				// If the element was set to 0 with this operation, we need to remove
 				// the pattern from the wave, and propagate the information
 				if (value[direction] == 0) {
-					add_to_propagator(y2, x2, *it);
-					wave.set(i2, *it, false);
+					add_to_propagator(y2, x2, pattern);
+					wave.set(i2, pattern, false);
 				}
 			}
 		}
diff --git a/modules/wfc/propagator.h b/modules/wfc/propagator.h
index 3f6e3dd72..bde774066 100644
--- a/modules/wfc/propagator.h
+++ b/modules/wfc/propagator.h
@@ -5,13 +5,13 @@
 #include "direction.h"
 #include <array>
 #include <tuple>
-#include <vector>
+#include "core/vector.h"
 
 class Wave;
 
 class Propagator {
 public:
-	using PropagatorState = std::vector<std::array<std::vector<unsigned>, 4>>;
+	using PropagatorState = Vector<std::array<Vector<unsigned>, 4>>;
 
 private:
 	const std::size_t patterns_size;
@@ -26,7 +26,7 @@ private:
 	// All the tuples (y, x, pattern) that should be propagated.
 	// The tuple should be propagated when wave.get(y, x, pattern) is set to
 	// false.
-	std::vector<std::tuple<unsigned, unsigned, unsigned>> propagating;
+	Vector<std::tuple<unsigned, unsigned, unsigned>> propagating;
 
 	// compatible.get(y, x, pattern)[direction] contains the number of patterns
 	// present in the wave that can be placed in the cell next to (y,x) in the
@@ -53,7 +53,8 @@ public:
 		// All the direction are set to 0, since the pattern cannot be set in (y,x).
 		std::array<int, 4> temp = {};
 		compatible.get(y, x, pattern) = temp;
-		propagating.emplace_back(y, x, pattern);
+		
+		propagating.push_back(std::tuple<unsigned, unsigned, unsigned>(y, x, pattern));
 	}
 
 	void propagate(Wave &wave);
diff --git a/modules/wfc/tiling_wfc.h b/modules/wfc/tiling_wfc.h
index 947d588ef..4585117a3 100644
--- a/modules/wfc/tiling_wfc.h
+++ b/modules/wfc/tiling_wfc.h
@@ -2,7 +2,7 @@
 #define FAST_WFC_TILING_WFC_HPP_
 
 #include <unordered_map>
-#include <vector>
+#include "core/vector.h"
 
 #include "array_2d.h"
 #include "wfc.h"
@@ -40,13 +40,13 @@ constexpr unsigned nb_of_possible_orientations(const Symmetry &symmetry) {
 // A tile that can be placed on the board.
 template <typename T>
 struct Tile {
-	std::vector<Array2D<T>> data; // The different orientations of the tile
+	Vector<Array2D<T>> data; // The different orientations of the tile
 	Symmetry symmetry; // The symmetry of the tile
 	double weight; // Its weight on the distribution of presence of tiles
 
 	// Generate the map associating an orientation id to the orientation
 	// id obtained when rotating 90° anticlockwise the tile.
-	static std::vector<unsigned> generate_rotation_map(const Symmetry &symmetry) {
+	static Vector<unsigned> generate_rotation_map(const Symmetry &symmetry) {
 		switch (symmetry) {
 			case Symmetry::X:
 				return { 0 };
@@ -64,7 +64,7 @@ struct Tile {
 
 	// Generate the map associating an orientation id to the orientation
 	// id obtained when reflecting the tile along the x axis.
-	static std::vector<unsigned> generate_reflection_map(const Symmetry &symmetry) {
+	static Vector<unsigned> generate_reflection_map(const Symmetry &symmetry) {
 		switch (symmetry) {
 			case Symmetry::X:
 				return { 0 };
@@ -87,12 +87,12 @@ struct Tile {
 	// Actions 0, 1, 2, and 3 are 0°, 90°, 180°, and 270° anticlockwise rotations.
 	// Actions 4, 5, 6, and 7 are actions 0, 1, 2, and 3 preceded by a reflection
 	// on the x axis.
-	static std::vector<std::vector<unsigned>> generate_action_map(const Symmetry &symmetry) {
-		std::vector<unsigned> rotation_map = generate_rotation_map(symmetry);
-		std::vector<unsigned> reflection_map = generate_reflection_map(symmetry);
+	static Vector<Vector<unsigned>> generate_action_map(const Symmetry &symmetry) {
+		Vector<unsigned> rotation_map = generate_rotation_map(symmetry);
+		Vector<unsigned> reflection_map = generate_reflection_map(symmetry);
 		size_t size = rotation_map.size();
-		std::vector<std::vector<unsigned>> action_map(8,
-				std::vector<unsigned>(size));
+		Vector<Vector<unsigned>> action_map(8,
+				Vector<unsigned>(size));
 		for (size_t i = 0; i < size; ++i) {
 			action_map[0][i] = i;
 		}
@@ -114,9 +114,9 @@ struct Tile {
 	}
 
 	// Generate all distincts rotations of a 2D array given its symmetries;
-	static std::vector<Array2D<T>> generate_oriented(Array2D<T> data,
+	static Vector<Array2D<T>> generate_oriented(Array2D<T> data,
 			Symmetry symmetry) {
-		std::vector<Array2D<T>> oriented;
+		Vector<Array2D<T>> oriented;
 		oriented.push_back(data);
 
 		switch (symmetry) {
@@ -148,7 +148,7 @@ struct Tile {
 
 	// Create a tile with its differents orientations, its symmetries and its
 	// weight on the distribution of tiles.
-	Tile(std::vector<Array2D<T>> data, Symmetry symmetry, double weight) :
+	Tile(Vector<Array2D<T>> data, Symmetry symmetry, double weight) :
 			data(data), symmetry(symmetry), weight(weight) {}
 
 	// Create a tile with its base orientation, its symmetries and its
@@ -166,9 +166,9 @@ struct TilingWFCOptions {
 template <typename T>
 class TilingWFC {
 private:
-	std::vector<Tile<T>> tiles;
-	std::vector<std::pair<unsigned, unsigned>> id_to_oriented_tile;
-	std::vector<std::vector<unsigned>> oriented_tile_ids;
+	Vector<Tile<T>> tiles;
+	Vector<std::pair<unsigned, unsigned>> id_to_oriented_tile;
+	Vector<Vector<unsigned>> oriented_tile_ids;
 
 	TilingWFCOptions options;
 
@@ -180,11 +180,11 @@ public:
 
 private:
 	// Generate mapping from id to oriented tiles and vice versa.
-	static std::pair<std::vector<std::pair<unsigned, unsigned>>,
-			std::vector<std::vector<unsigned>>>
-	generate_oriented_tile_ids(const std::vector<Tile<T>> &tiles) {
-		std::vector<std::pair<unsigned, unsigned>> id_to_oriented_tile;
-		std::vector<std::vector<unsigned>> oriented_tile_ids;
+	static std::pair<Vector<std::pair<unsigned, unsigned>>,
+			Vector<Vector<unsigned>>>
+	generate_oriented_tile_ids(const Vector<Tile<T>> &tiles) {
+		Vector<std::pair<unsigned, unsigned>> id_to_oriented_tile;
+		Vector<Vector<unsigned>> oriented_tile_ids;
 
 		unsigned id = 0;
 		for (unsigned i = 0; i < tiles.size(); i++) {
@@ -200,38 +200,34 @@ private:
 	}
 
 	// Generate the propagator which will be used in the wfc algorithm.
-	static std::vector<std::array<std::vector<unsigned>, 4>> generate_propagator(
-			const std::vector<std::tuple<unsigned, unsigned, unsigned, unsigned>>
-					&neighbors,
-			std::vector<Tile<T>> tiles,
-			std::vector<std::pair<unsigned, unsigned>> id_to_oriented_tile,
-			std::vector<std::vector<unsigned>> oriented_tile_ids) {
+	static Vector<std::array<Vector<unsigned>, 4>> generate_propagator(
+			const Vector<std::tuple<unsigned, unsigned, unsigned, unsigned>> &neighbors,
+			Vector<Tile<T>> tiles,
+			Vector<std::pair<unsigned, unsigned>> id_to_oriented_tile,
+			Vector<Vector<unsigned>> oriented_tile_ids) {
+
 		size_t nb_oriented_tiles = id_to_oriented_tile.size();
-		std::vector<std::array<std::vector<bool>, 4>> dense_propagator(
-				nb_oriented_tiles, { std::vector<bool>(nb_oriented_tiles, false), std::vector<bool>(nb_oriented_tiles, false), std::vector<bool>(nb_oriented_tiles, false), std::vector<bool>(nb_oriented_tiles, false) });
+		Vector<std::array<Vector<bool>, 4>> dense_propagator(
+				nb_oriented_tiles, { Vector<bool>(nb_oriented_tiles, false), Vector<bool>(nb_oriented_tiles, false), Vector<bool>(nb_oriented_tiles, false), Vector<bool>(nb_oriented_tiles, false) });
 
 		for (auto neighbor : neighbors) {
 			unsigned tile1 = std::get<0>(neighbor);
 			unsigned orientation1 = std::get<1>(neighbor);
 			unsigned tile2 = std::get<2>(neighbor);
 			unsigned orientation2 = std::get<3>(neighbor);
-			std::vector<std::vector<unsigned>> action_map1 =
+			Vector<Vector<unsigned>> action_map1 =
 					Tile<T>::generate_action_map(tiles[tile1].symmetry);
-			std::vector<std::vector<unsigned>> action_map2 =
+			Vector<Vector<unsigned>> action_map2 =
 					Tile<T>::generate_action_map(tiles[tile2].symmetry);
 
 			auto add = [&](unsigned action, unsigned direction) {
 				unsigned temp_orientation1 = action_map1[action][orientation1];
 				unsigned temp_orientation2 = action_map2[action][orientation2];
-				unsigned oriented_tile_id1 =
-						oriented_tile_ids[tile1][temp_orientation1];
-				unsigned oriented_tile_id2 =
-						oriented_tile_ids[tile2][temp_orientation2];
-				dense_propagator[oriented_tile_id1][direction][oriented_tile_id2] =
-						true;
+				unsigned oriented_tile_id1 = oriented_tile_ids[tile1][temp_orientation1];
+				unsigned oriented_tile_id2 = oriented_tile_ids[tile2][temp_orientation2];
+				dense_propagator[oriented_tile_id1][direction][oriented_tile_id2] = true;
 				direction = get_opposite_direction(direction);
-				dense_propagator[oriented_tile_id2][direction][oriented_tile_id1] =
-						true;
+				dense_propagator[oriented_tile_id2][direction][oriented_tile_id1] = true;
 			};
 
 			add(0, 2);
@@ -244,8 +240,8 @@ private:
 			add(7, 0);
 		}
 
-		std::vector<std::array<std::vector<unsigned>, 4>> propagator(
-				nb_oriented_tiles);
+		Vector<std::array<Vector<unsigned>, 4>> propagator(nb_oriented_tiles);
+
 		for (size_t i = 0; i < nb_oriented_tiles; ++i) {
 			for (size_t j = 0; j < nb_oriented_tiles; ++j) {
 				for (size_t d = 0; d < 4; ++d) {
@@ -260,14 +256,15 @@ private:
 	}
 
 	// Get probability of presence of tiles.
-	static std::vector<double>
-	get_tiles_weights(const std::vector<Tile<T>> &tiles) {
-		std::vector<double> frequencies;
+	static Vector<double> get_tiles_weights(const Vector<Tile<T>> &tiles) {
+		Vector<double> frequencies;
+
 		for (size_t i = 0; i < tiles.size(); ++i) {
 			for (size_t j = 0; j < tiles[i].data.size(); ++j) {
 				frequencies.push_back(tiles[i].weight / tiles[i].data.size());
 			}
 		}
+
 		return frequencies;
 	}
 
@@ -275,18 +272,19 @@ private:
 	Array2D<T> id_to_tiling(Array2D<unsigned> ids) {
 		unsigned size = tiles[0].data[0].height;
 		Array2D<T> tiling(size * ids.height, size * ids.width);
+
 		for (unsigned i = 0; i < ids.height; i++) {
 			for (unsigned j = 0; j < ids.width; j++) {
-				std::pair<unsigned, unsigned> oriented_tile =
-						id_to_oriented_tile[ids.get(i, j)];
+				std::pair<unsigned, unsigned> oriented_tile = id_to_oriented_tile[ids.get(i, j)];
+
 				for (unsigned y = 0; y < size; y++) {
 					for (unsigned x = 0; x < size; x++) {
-						tiling.get(i * size + y, j * size + x) =
-								tiles[oriented_tile.first].data[oriented_tile.second].get(y, x);
+						tiling.get(i * size + y, j * size + x) = tiles[oriented_tile.first].data[oriented_tile.second].get(y, x);
 					}
 				}
 			}
 		}
+
 		return tiling;
 	}
 
@@ -301,8 +299,8 @@ private:
 public:
 	// Construct the TilingWFC class to generate a tiled image.
 	TilingWFC(
-			const std::vector<Tile<T>> &tiles,
-			const std::vector<std::tuple<unsigned, unsigned, unsigned, unsigned>>
+			const Vector<Tile<T>> &tiles,
+			const Vector<std::tuple<unsigned, unsigned, unsigned, unsigned>>
 					&neighbors,
 			const unsigned height, const unsigned width,
 			const TilingWFCOptions &options, int seed) :
diff --git a/modules/wfc/wave.cpp b/modules/wfc/wave.cpp
index 390e97432..b82dc4499 100644
--- a/modules/wfc/wave.cpp
+++ b/modules/wfc/wave.cpp
@@ -5,27 +5,31 @@
 namespace {
 
 // Return distribution * log(distribution).
-std::vector<double> get_plogp(const std::vector<double> &distribution) {
-	std::vector<double> plogp;
+Vector<double> get_plogp(const Vector<double> &distribution) {
+	Vector<double> plogp;
+
 	for (unsigned i = 0; i < distribution.size(); i++) {
 		plogp.push_back(distribution[i] * log(distribution[i]));
 	}
+
 	return plogp;
 }
 
 // Return min(v) / 2.
-double get_min_abs_half(const std::vector<double> &v) {
+double get_min_abs_half(const Vector<double> &v) {
 	double min_abs_half = std::numeric_limits<double>::infinity();
+
 	for (unsigned i = 0; i < v.size(); i++) {
 		min_abs_half = std::min(min_abs_half, std::abs(v[i] / 2.0));
 	}
+
 	return min_abs_half;
 }
 
 } // namespace
 
 Wave::Wave(unsigned height, unsigned width,
-		const std::vector<double> &patterns_frequencies) :
+		const Vector<double> &patterns_frequencies) :
 		patterns_frequencies(patterns_frequencies),
 		plogp_patterns_frequencies(get_plogp(patterns_frequencies)),
 		min_abs_half_plogp(get_min_abs_half(plogp_patterns_frequencies)),
@@ -35,38 +39,52 @@ Wave::Wave(unsigned height, unsigned width,
 		width(width),
 		height(height),
 		size(height * width) {
+
 	// Initialize the memoisation of entropy.
 	double base_entropy = 0;
 	double base_s = 0;
+
 	for (unsigned i = 0; i < nb_patterns; i++) {
 		base_entropy += plogp_patterns_frequencies[i];
 		base_s += patterns_frequencies[i];
 	}
+
 	double log_base_s = log(base_s);
 	double entropy_base = log_base_s - base_entropy / base_s;
-	memoisation.plogp_sum = std::vector<double>(width * height, base_entropy);
-	memoisation.sum = std::vector<double>(width * height, base_s);
-	memoisation.log_sum = std::vector<double>(width * height, log_base_s);
-	memoisation.nb_patterns =
-			std::vector<unsigned>(width * height, static_cast<unsigned>(nb_patterns));
-	memoisation.entropy = std::vector<double>(width * height, entropy_base);
+
+	memoisation.plogp_sum.resize(width * height);
+	memoisation.plogp_sum.fill(base_entropy);
+
+	memoisation.sum.resize(width * height);
+	memoisation.sum.fill(base_s);
+
+	memoisation.log_sum.resize(width * height);
+	memoisation.log_sum.fill(log_base_s);
+
+	memoisation.nb_patterns.resize(width * height);
+	memoisation.nb_patterns.fill(static_cast<unsigned>(nb_patterns));
+
+	memoisation.entropy.resize(width * height);
+	memoisation.entropy.fill(entropy_base);
 }
 
 void Wave::set(unsigned index, unsigned pattern, bool value) {
 	bool old_value = data.get(index, pattern);
+
 	// If the value isn't changed, nothing needs to be done.
 	if (old_value == value) {
 		return;
 	}
+
 	// Otherwise, the memoisation should be updated.
 	data.get(index, pattern) = value;
-	memoisation.plogp_sum[index] -= plogp_patterns_frequencies[pattern];
-	memoisation.sum[index] -= patterns_frequencies[pattern];
-	memoisation.log_sum[index] = log(memoisation.sum[index]);
-	memoisation.nb_patterns[index]--;
-	memoisation.entropy[index] = memoisation.log_sum[index] - memoisation.plogp_sum[index] / memoisation.sum[index];
-	// If there is no patterns possible in the cell, then there is a
-	// contradiction.
+	memoisation.plogp_sum.write[index] -= plogp_patterns_frequencies[pattern];
+	memoisation.sum.write[index] -= patterns_frequencies[pattern];
+	memoisation.log_sum.write[index] = log(memoisation.sum[index]);
+	memoisation.nb_patterns.write[index]--;
+	memoisation.entropy.write[index] = memoisation.log_sum[index] - memoisation.plogp_sum[index] / memoisation.sum[index];
+
+	// If there is no patterns possible in the cell, then there is a contradiction.
 	if (memoisation.nb_patterns[index] == 0) {
 		is_impossible = true;
 	}
@@ -87,6 +105,7 @@ int Wave::get_min_entropy(std::minstd_rand &gen) const {
 		// If the cell is decided, we do not compute the entropy (which is equal
 		// to 0).
 		double nb_patterns_local = memoisation.nb_patterns[i];
+		
 		if (nb_patterns_local == 1) {
 			continue;
 		}
diff --git a/modules/wfc/wave.h b/modules/wfc/wave.h
index 4ad71cf2e..2fb23b032 100644
--- a/modules/wfc/wave.h
+++ b/modules/wfc/wave.h
@@ -3,18 +3,18 @@
 
 #include "array_2d.h"
 #include <random>
-#include <vector>
+#include "core/vector.h"
 
 // Struct containing the values needed to compute the entropy of all the cells.
 // This struct is updated every time the wave is changed.
 // p'(pattern) is equal to patterns_frequencies[pattern] if wave.get(cell,
 // pattern) is set to true, otherwise 0.
 struct EntropyMemoisation {
-	std::vector<double> plogp_sum; // The sum of p'(pattern)// log(p'(pattern)).
-	std::vector<double> sum; // The sum of p'(pattern).
-	std::vector<double> log_sum; // The log of sum.
-	std::vector<unsigned> nb_patterns; // The number of patterns present
-	std::vector<double> entropy; // The entropy of the cell.
+	Vector<double> plogp_sum; // The sum of p'(pattern)// log(p'(pattern)).
+	Vector<double> sum; // The sum of p'(pattern).
+	Vector<double> log_sum; // The log of sum.
+	Vector<unsigned> nb_patterns; // The number of patterns present
+	Vector<double> entropy; // The entropy of the cell.
 };
 
 // Contains the pattern possibilities in every cell.
@@ -22,10 +22,10 @@ struct EntropyMemoisation {
 class Wave {
 private:
 	// The patterns frequencies p given to wfc.
-	const std::vector<double> patterns_frequencies;
+	const Vector<double> patterns_frequencies;
 
 	// The precomputation of p * log(p).
-	const std::vector<double> plogp_patterns_frequencies;
+	const Vector<double> plogp_patterns_frequencies;
 
 	// The precomputation of min (p * log(p)) / 2.
 	// This is used to define the maximum value of the noise.
@@ -50,7 +50,7 @@ public:
 	const unsigned size;
 
 	// Initialize the wave with every cell being able to have every pattern.
-	Wave(unsigned height, unsigned width, const std::vector<double> &patterns_frequencies);
+	Wave(unsigned height, unsigned width, const Vector<double> &patterns_frequencies);
 
 	// Return true if pattern can be placed in cell index.
 	bool get(unsigned index, unsigned pattern) const {
diff --git a/modules/wfc/wfc.cpp b/modules/wfc/wfc.cpp
index 4894c7c14..fd73ed7cf 100644
--- a/modules/wfc/wfc.cpp
+++ b/modules/wfc/wfc.cpp
@@ -3,15 +3,18 @@
 
 namespace {
 // Normalize a vector so the sum of its elements is equal to 1.0f
-std::vector<double> &normalize(std::vector<double> &v) {
+Vector<double> &normalize(Vector<double> &v) {
 	double sum_weights = 0.0;
-	for (double weight : v) {
-		sum_weights += weight;
+	int size = v.size();
+	const double* vpr = v.ptr();
+	for (int i = 0; i < size; ++i) {
+		sum_weights += vpr[i];
 	}
 
+	double* vpw = v.ptrw();
 	double inv_sum_weights = 1.0 / sum_weights;
-	for (double &weight : v) {
-		weight *= inv_sum_weights;
+	for (int i = 0; i < size; ++i) {
+		vpw[i] *= inv_sum_weights;
 	}
 
 	return v;
@@ -20,10 +23,11 @@ std::vector<double> &normalize(std::vector<double> &v) {
 
 Array2D<unsigned> WFC::wave_to_output() const {
 	Array2D<unsigned> output_patterns(wave.height, wave.width);
+
 	for (unsigned i = 0; i < wave.size; i++) {
 		for (unsigned k = 0; k < nb_patterns; k++) {
 			if (wave.get(i, k)) {
-				output_patterns.data[i] = k;
+				output_patterns.data.write[i] = k;
 			}
 		}
 	}
@@ -31,7 +35,7 @@ Array2D<unsigned> WFC::wave_to_output() const {
 }
 
 WFC::WFC(bool periodic_output, int seed,
-		std::vector<double> patterns_frequencies,
+		Vector<double> patterns_frequencies,
 		Propagator::PropagatorState propagator, unsigned wave_height,
 		unsigned wave_width) :
 		gen(seed), patterns_frequencies(normalize(patterns_frequencies)), wave(wave_height, wave_width, patterns_frequencies), nb_patterns(propagator.size()), propagator(wave.height, wave.width, periodic_output, propagator) {}
diff --git a/modules/wfc/wfc.h b/modules/wfc/wfc.h
index 0914ec9d0..25aa9fb45 100644
--- a/modules/wfc/wfc.h
+++ b/modules/wfc/wfc.h
@@ -15,7 +15,7 @@ private:
 	std::minstd_rand gen;
 
 	// The distribution of the patterns as given in input.
-	const std::vector<double> patterns_frequencies;
+	const Vector<double> patterns_frequencies;
 
 	Wave wave;
 
@@ -32,7 +32,7 @@ private:
 
 public:
 	// Basic constructor initializing the algorithm.
-	WFC(bool periodic_output, int seed, std::vector<double> patterns_frequencies,
+	WFC(bool periodic_output, int seed, Vector<double> patterns_frequencies,
 			Propagator::PropagatorState propagator, unsigned wave_height,
 			unsigned wave_width);