diff --git a/engibench/problems/thermoelastic2d/model/fea_model.py b/engibench/problems/thermoelastic2d/model/fea_model.py
index af4eb673..ab84490a 100644
--- a/engibench/problems/thermoelastic2d/model/fea_model.py
+++ b/engibench/problems/thermoelastic2d/model/fea_model.py
@@ -119,7 +119,7 @@ def run(self, bcs: dict[str, Any], x_init: np.ndarray | None = None) -> dict[str
         Args:
             bcs (dict[str, any]): A dictionary containing boundary conditions and problem parameters.
                 Expected keys include:
-                    - 'volfrac' (float): Target volume fraction.
+                    - 'volume_fraction_target' (float): Target volume fraction.
                     - 'fixed_elements' (np.ndarray): NxN binary array encoding the location of fixed elements.
                     - 'force_elements_x' (np.ndarray): NxN binary array encoding the location of loaded elements in the x direction.
                     - 'force_elements_y' (np.ndarray): NxN binary array encoding the location of loaded elements in the y direction.
@@ -146,7 +146,7 @@ def run(self, bcs: dict[str, Any], x_init: np.ndarray | None = None) -> dict[str
         nelx = fe_h - 1
         nely = fe_w - 1
 
-        volfrac = bcs["volfrac"]
+        volfrac = bcs["volume_fraction_target"]
         n = nely * nelx  # Total number of elements
 
         # OptiSteps records
@@ -285,7 +285,7 @@ def run(self, bcs: dict[str, Any], x_init: np.ndarray | None = None) -> dict[str
                 return {
                     "structural_compliance": f0valm,
                     "thermal_compliance": f0valt,
-                    "volume_fraction": vf_error,
+                    "volume_fraction_error": vf_error,
                 }
             vf_error = np.abs(np.mean(x) - volfrac)
             obj_values = np.array([f0valm, f0valt, vf_error])
@@ -353,7 +353,7 @@ def run(self, bcs: dict[str, Any], x_init: np.ndarray | None = None) -> dict[str
             "bcs": bcs,
             "structural_compliance": f0valm,
             "thermal_compliance": f0valt,
-            "volume_fraction": vf_error,
+            "volume_fraction_error": vf_error,
             "opti_steps": opti_steps,
         }
 
@@ -372,7 +372,7 @@ def run(self, bcs: dict[str, Any], x_init: np.ndarray | None = None) -> dict[str
         "fixed_elements": [lci[21], lci[32], lci[43]],
         "force_elements_y": [bri[31]],
         "heatsink_elements": [lci[31], lci[32], lci[33]],
-        "volfrac": 0.2,
+        "volume_fraction_target": 0.2,
         "rmin": 1.1,
         "weight": 1.0,  # 1.0 for pure structural, 0.0 for pure thermal
     }
diff --git a/engibench/problems/thermoelastic2d/v0.py b/engibench/problems/thermoelastic2d/v0.py
index 00e34d9b..61c701d8 100644
--- a/engibench/problems/thermoelastic2d/v0.py
+++ b/engibench/problems/thermoelastic2d/v0.py
@@ -42,7 +42,7 @@ class ThermoElastic2D(Problem[npt.NDArray]):
     objectives: tuple[tuple[str, ObjectiveDirection], ...] = (
         ("structural_compliance", ObjectiveDirection.MINIMIZE),
         ("thermal_compliance", ObjectiveDirection.MINIMIZE),
-        ("volume_fraction", ObjectiveDirection.MINIMIZE),
+        ("volume_fraction_error", ObjectiveDirection.MINIMIZE),
     )
 
     @dataclass
@@ -65,7 +65,7 @@ class Conditions:
             default_factory=lambda: HEATSINK_ELEMENTS
         )
         """Binary NxN matrix specifying elements that have a heat sink"""
-        volfrac: Annotated[float, bounded(lower=0.0, upper=1.0).category(THEORY)] = 0.3
+        volume_fraction_target: Annotated[float, bounded(lower=0.0, upper=1.0).category(THEORY)] = 0.3
         """Target volume fraction for the volume fraction constraint"""
         rmin: Annotated[
             float, bounded(lower=1.0).category(THEORY), bounded(lower=0.0, upper=3.0).warning().category(IMPL)
@@ -129,7 +129,7 @@ def simulate(self, design: npt.NDArray, config: dict[str, Any] | None = None) ->
                     boundary_dict[key] = value
 
         results = FeaModel(plot=False, eval_only=True).run(boundary_dict, x_init=design)
-        return np.array([results["structural_compliance"], results["thermal_compliance"], results["volume_fraction"]])
+        return np.array([results["structural_compliance"], results["thermal_compliance"], results["volume_fraction_error"]])
 
     def optimize(
         self, starting_point: npt.NDArray, config: dict[str, Any] | None = None