diff --git a/projection_utils.py b/projection_utils.py
index 7b5d7749ad93382c25fe1e00ff7c43839f29f37c..1bf0eedc7aec1b1d06cb7aeefbde6e8daff3cacb 100644
--- a/projection_utils.py
+++ b/projection_utils.py
@@ -4,6 +4,7 @@
 
 """
 
+from functools import cache
 from typing import Tuple
 import numpy as np
 from numpy import ndarray
@@ -16,6 +17,78 @@ from Initial_Condition import InitialCondition
 x = Symbol('x')
 
 
+class Mesh:
+    """Class for mesh/grid.
+
+    Each cell is characterized by its center.
+
+    Attributes
+    ----------
+    num_grid_cells : int
+        Number of cells in the mesh (ghost cells notwithstanding). Usually
+        exponential of 2.
+    bounds : Tuple[float, float]
+        Left and right boundary of the mesh interval.
+    interval_len : float
+        Length of the mesh interval.
+    cell_len : float
+        Length of a mesh cell.
+    cells : ndarray
+        Array of cell centers in mesh.
+
+    """
+
+    def __init__(self, num_grid_cells: int, num_ghost_cells: int,
+                 left_bound: float, right_bound: float) -> None:
+        """Initialize Mesh.
+
+        Parameters
+        ----------
+        num_grid_cells : int
+            Number of cells in the mesh (ghost cells notwithstanding). Usually
+            exponential of 2.
+        num_ghost_cells : int
+            Number of ghost cells on each side of the mesh.
+        left_bound : float
+            Left boundary of the mesh interval.
+        right_bound : float
+            Right boundary of the mesh interval.
+        """
+        self._num_grid_cells = num_grid_cells
+        self._num_ghost_cells = num_ghost_cells
+        self._left_bound = left_bound
+        self._right_bound = right_bound
+
+    @property
+    def num_grid_cells(self) -> int:
+        """Return number of grid cells."""
+        return self._num_grid_cells
+
+    @property
+    def bounds(self) -> Tuple[float, float]:
+        """Return left and right boundary of the mesh interval."""
+        return self._left_bound, self._right_bound
+
+    @property
+    def interval_len(self) -> float:
+        """Return the length of the mesh interval."""
+        return self._right_bound - self._left_bound
+
+    @property
+    def cell_len(self) -> float:
+        """Return the length of a mesh cell."""
+        return self.interval_len/self.num_grid_cells
+
+    @property
+    @cache
+    def cells(self) -> ndarray:
+        """Return the cell centers of the mesh (including ghost cells)."""
+        return np.arange(
+            self._left_bound - (self._num_ghost_cells*2-1)/2*self.cell_len,
+            self._right_bound + (self._num_ghost_cells*2+1)/2*self.cell_len,
+            self.cell_len)
+
+
 def calculate_approximate_solution(
         projection: ndarray, points: ndarray, polynomial_degree: int,
         basis: ndarray) -> ndarray: