diff --git a/src/mechaphlowers/api/section_study.py b/src/mechaphlowers/api/section_study.py
index 6d920540..f461ddb6 100644
--- a/src/mechaphlowers/api/section_study.py
+++ b/src/mechaphlowers/api/section_study.py
@@ -32,6 +32,7 @@
 
 if TYPE_CHECKING:
     from mechaphlowers.core.models.cable.thermal import ThermalEngine
+    from mechaphlowers.core.models.estimation import EstimationEngine
     from mechaphlowers.core.models.guying import Guying
     from mechaphlowers.plotting.plot import PlotEngine
 from mechaphlowers.utils import arr
@@ -94,6 +95,7 @@ def __init__(
         self._plot_engine: PlotEngine | None = None
         self._thermal_engine: ThermalEngine | None = None
         self._guying: Guying | None = None
+        self._estimation_engine: EstimationEngine | None = None
         self._intermediate_memento: BalanceEngineMemento | None = None
 
     # ── Sub-engine properties ─────────────────────────────────────────────
@@ -132,6 +134,25 @@ def guying(self) -> Guying:
             self._guying = _G(self._balance_engine)
         return self._guying
 
+    @property
+    def estimation_engine(self) -> EstimationEngine:
+        """Lazy-loaded inverse estimation engine.
+
+        Returns an [`EstimationEngine`][mechaphlowers.core.models.estimation.EstimationEngine]
+        bound to this study, using Brent's method by default.
+        """
+        if self._estimation_engine is None:
+            from mechaphlowers.core.models.estimation import (
+                EstimationEngine as _EE,
+            )
+
+            self._estimation_engine = _EE(self)
+        return self._estimation_engine
+
+    @estimation_engine.setter
+    def estimation_engine(self, estimation_engine: EstimationEngine):
+        self._estimation_engine = estimation_engine
+
     @property
     def intermediate_memento(self) -> BalanceEngineMemento | None:
         """The memento captured after the intermediate warm-start solve, if any."""
diff --git a/src/mechaphlowers/core/models/estimation/__init__.py b/src/mechaphlowers/core/models/estimation/__init__.py
new file mode 100644
index 00000000..a89e4e25
--- /dev/null
+++ b/src/mechaphlowers/core/models/estimation/__init__.py
@@ -0,0 +1,23 @@
+# Copyright (c) 2026, RTE (http://www.rte-france.com)
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+# SPDX-License-Identifier: MPL-2.0
+
+from mechaphlowers.core.models.estimation.engine import EstimationEngine
+from mechaphlowers.core.models.estimation.methods import (
+    BisectionMethod,
+    BrentMethod,
+    NewtonMethod,
+    OptimizationMethod,
+)
+from mechaphlowers.core.models.estimation.result import EstimationResult
+
+__all__ = [
+    "EstimationEngine",
+    "EstimationResult",
+    "OptimizationMethod",
+    "BisectionMethod",
+    "BrentMethod",
+    "NewtonMethod",
+]
diff --git a/src/mechaphlowers/core/models/estimation/engine.py b/src/mechaphlowers/core/models/estimation/engine.py
new file mode 100644
index 00000000..f44fdd85
--- /dev/null
+++ b/src/mechaphlowers/core/models/estimation/engine.py
@@ -0,0 +1,285 @@
+# Copyright (c) 2026, RTE (http://www.rte-france.com)
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+# SPDX-License-Identifier: MPL-2.0
+
+from __future__ import annotations
+
+import logging
+from typing import TYPE_CHECKING, Callable
+
+import numpy as np
+
+from mechaphlowers.core.models.estimation.methods import (
+    BrentMethod,
+    OptimizationMethod,
+)
+from mechaphlowers.core.models.estimation.result import EstimationResult
+
+if TYPE_CHECKING:
+    from mechaphlowers.api.section_study import SectionStudy
+
+logger = logging.getLogger(__name__)
+
+
+class EstimationEngine:
+    """Generic inverse-problem solver built on top of SectionStudy.
+
+    Wraps a `SectionStudy` and an `OptimizationMethod` to find the value of
+    a physical variable (temperature, wind, load) that produces a target
+    distance to an obstacle.
+
+    The engine saves/restores the balance-engine state around every objective
+    evaluation so that the study is left unchanged after estimation.
+
+    Args:
+        study: The `SectionStudy` instance (must have been solved via
+            `solve_adjustment` beforehand).
+        method: The optimization algorithm to use. Defaults to `BrentMethod`.
+
+    Examples:
+        >>> engine = EstimationEngine(study, method=BrentMethod(tol=0.01))
+        >>> result = engine.estimate_temperature(
+        ...     span_index=0,
+        ...     obstacle_point=np.array([150.0, 0.0, 5.0]),
+        ...     target_distance=8.0,
+        ...     bounds=(0.0, 200.0),
+        ... )
+        >>> print(result.value, result.converged)
+    """
+
+    def __init__(
+        self,
+        study: SectionStudy,
+        method: OptimizationMethod | None = None,
+    ) -> None:
+        self._study = study
+        self._method: OptimizationMethod = method or BrentMethod()
+
+    @property
+    def method(self) -> OptimizationMethod:
+        return self._method
+
+    @method.setter
+    def method(self, value: OptimizationMethod) -> None:
+        self._method = value
+
+    def estimate(
+        self,
+        objective: Callable[[float], float],
+        bounds: tuple[float, float],
+    ) -> EstimationResult:
+        """Run the optimization on a generic objective function.
+
+        The objective must be a callable ``f(x) -> float`` where the root
+        ``f(x) = 0`` corresponds to the desired solution. State management
+        (save/restore) is the caller's responsibility when using this method
+        directly.
+
+        Args:
+            objective: Function to zero. Signature: ``(x: float) -> float``.
+            bounds: ``(lower, upper)`` search interval.
+
+        Returns:
+            EstimationResult with the solution.
+        """
+        return self._method.solve(objective, bounds)
+
+    def estimate_temperature(
+        self,
+        span_index: int,
+        obstacle_point: np.ndarray,
+        target_distance: float,
+        bounds: tuple[float, float] = (0.0, 200.0),
+        wind_pressure: float | None = None,
+        ice_thickness: float | None = None,
+    ) -> EstimationResult:
+        """Find the cable temperature that yields a target distance to an obstacle.
+
+        Args:
+            span_index: Index of the span where the obstacle is located.
+            obstacle_point: 3D coordinates of the obstacle point (shape ``(3,)``).
+            target_distance: Desired distance between cable and obstacle (meters).
+            bounds: Search interval for temperature in °C.
+            wind_pressure: Fixed wind pressure in Pa (optional).
+            ice_thickness: Fixed ice thickness in m (optional).
+
+        Returns:
+            EstimationResult with the temperature value.
+        """
+
+        def objective(temperature: float) -> float:
+            return self._distance_difference(
+                span_index=span_index,
+                obstacle_point=obstacle_point,
+                target_distance=target_distance,
+                new_temperature=temperature,
+                wind_pressure=wind_pressure,
+                ice_thickness=ice_thickness,
+            )
+
+        logger.info(
+            "Estimating temperature for target distance %.3f m on span %d",
+            target_distance,
+            span_index,
+        )
+        return self._method.solve(objective, bounds)
+
+    def estimate_wind(
+        self,
+        span_index: int,
+        obstacle_point: np.ndarray,
+        target_distance: float,
+        bounds: tuple[float, float] = (0.0, 2000.0),
+        new_temperature: float | None = None,
+        ice_thickness: float | None = None,
+    ) -> EstimationResult:
+        """Find the wind pressure that yields a target distance to an obstacle.
+
+        Args:
+            span_index: Index of the span where the obstacle is located.
+            obstacle_point: 3D coordinates of the obstacle point (shape ``(3,)``).
+            target_distance: Desired distance between cable and obstacle (meters).
+            bounds: Search interval for wind pressure in Pa.
+            new_temperature: Fixed temperature in °C (optional).
+            ice_thickness: Fixed ice thickness in m (optional).
+
+        Returns:
+            EstimationResult with the wind pressure value.
+        """
+
+        def objective(wind: float) -> float:
+            return self._distance_difference(
+                span_index=span_index,
+                obstacle_point=obstacle_point,
+                target_distance=target_distance,
+                wind_pressure=wind,
+                new_temperature=new_temperature,
+                ice_thickness=ice_thickness,
+            )
+
+        logger.info(
+            "Estimating wind pressure for target distance %.3f m on span %d",
+            target_distance,
+            span_index,
+        )
+        return self._method.solve(objective, bounds)
+
+    def estimate_load(
+        self,
+        span_index: int,
+        obstacle_point: np.ndarray,
+        target_distance: float,
+        load_position_distance: float,
+        bounds: tuple[float, float] = (0.0, 100.0),
+        new_temperature: float | None = None,
+        wind_pressure: float | None = None,
+        ice_thickness: float | None = None,
+    ) -> EstimationResult:
+        """Find the load mass that yields a target distance to an obstacle.
+
+        Args:
+            span_index: Index of the span where the obstacle is located.
+            obstacle_point: 3D coordinates of the obstacle point (shape ``(3,)``).
+            target_distance: Desired distance between cable and obstacle (meters).
+            load_position_distance: Position of the load along the span (meters).
+            bounds: Search interval for load mass in kg.
+            new_temperature: Fixed temperature in °C (optional).
+            wind_pressure: Fixed wind pressure in Pa (optional).
+            ice_thickness: Fixed ice thickness in m (optional).
+
+        Returns:
+            EstimationResult with the load mass value.
+        """
+
+        def objective(load_mass: float) -> float:
+            return self._distance_difference_with_load(
+                span_index=span_index,
+                obstacle_point=obstacle_point,
+                target_distance=target_distance,
+                load_position_distance=load_position_distance,
+                load_mass=load_mass,
+                new_temperature=new_temperature,
+                wind_pressure=wind_pressure,
+                ice_thickness=ice_thickness,
+            )
+
+        logger.info(
+            "Estimating load mass for target distance %.3f m on span %d",
+            target_distance,
+            span_index,
+        )
+        return self._method.solve(objective, bounds)
+
+    # ── Private helpers ───────────────────────────────────────────────────
+
+    def _distance_difference(
+        self,
+        span_index: int,
+        obstacle_point: np.ndarray,
+        target_distance: float,
+        wind_pressure: float | None = None,
+        ice_thickness: float | None = None,
+        new_temperature: float | None = None,
+    ) -> float:
+        """Compute ``distance(x) - target`` with state save/restore.
+
+        Solves change-state with the given parameters, computes the distance
+        to the obstacle, then restores the engine to its original state.
+        """
+        memento = self._study.save_state()
+        try:
+            self._study.solve_change_state(
+                wind_pressure=wind_pressure,
+                ice_thickness=ice_thickness,
+                new_temperature=new_temperature,
+            )
+            distance_result = self._study.position_engine.point_distance(
+                span_index, obstacle_point
+            )
+            distance = distance_result.distance_3d
+        finally:
+            self._study.restore_state(memento)
+
+        return distance - target_distance
+
+    def _distance_difference_with_load(
+        self,
+        span_index: int,
+        obstacle_point: np.ndarray,
+        target_distance: float,
+        load_position_distance: float,
+        load_mass: float,
+        new_temperature: float | None = None,
+        wind_pressure: float | None = None,
+        ice_thickness: float | None = None,
+    ) -> float:
+        """Compute distance difference between target distance and distance with a given load mass.
+
+        Saves state before computiong and restores it afterwards."""
+        memento = self._study.save_state()
+        try:
+            # Build load arrays for single point load
+            n_spans = len(
+                self._study.balance_engine.section_array.data.span_length
+            )
+            load_positions = np.zeros(n_spans)
+            load_masses = np.zeros(n_spans)
+            load_positions[span_index] = load_position_distance
+            load_masses[span_index] = load_mass
+
+            self._study.set_loads(load_positions, load_masses)
+            self._study.solve_change_state(
+                wind_pressure=wind_pressure,
+                ice_thickness=ice_thickness,
+                new_temperature=new_temperature,
+            )
+            distance_result = self._study.position_engine.point_distance(
+                span_index, obstacle_point
+            )
+            distance = distance_result.distance_3d
+        finally:
+            self._study.restore_state(memento)
+
+        return distance - target_distance
diff --git a/src/mechaphlowers/core/models/estimation/methods.py b/src/mechaphlowers/core/models/estimation/methods.py
new file mode 100644
index 00000000..348f347d
--- /dev/null
+++ b/src/mechaphlowers/core/models/estimation/methods.py
@@ -0,0 +1,220 @@
+# Copyright (c) 2026, RTE (http://www.rte-france.com)
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+# SPDX-License-Identifier: MPL-2.0
+
+from __future__ import annotations
+
+import logging
+from typing import Callable, Protocol
+
+from mechaphlowers.core.models.estimation.result import EstimationResult
+
+logger = logging.getLogger(__name__)
+
+
+class OptimizationMethod(Protocol):
+    """Protocol for root-finding / optimization algorithms.
+
+    Implementations solve ``objective(x) = 0`` within given bounds.
+    """
+
+    def solve(
+        self,
+        objective: Callable[[float], float],
+        bounds: tuple[float, float],
+    ) -> EstimationResult: ...
+
+
+class BisectionMethod:
+    """Simple bisection root-finding.
+
+    Requires that ``objective`` changes sign over ``bounds``.
+
+    Args:
+        tol: Absolute tolerance on the root. Defaults to 1e-3.
+        maxiter: Maximum number of iterations. Defaults to 50.
+    """
+
+    def __init__(self, tol: float = 1e-3, maxiter: int = 50) -> None:
+        self.tol = tol
+        self.maxiter = maxiter
+
+    def solve(
+        self,
+        objective: Callable[[float], float],
+        bounds: tuple[float, float],
+    ) -> EstimationResult:
+        a, b = bounds
+        fa = objective(a)
+        fb = objective(b)
+        iterations = 2
+
+        if fa * fb > 0:
+            logger.warning(
+                "Bisection: objective does not change sign over bounds "
+                f"[{a}, {b}] (f(a)={fa:.3e}, f(b)={fb:.3e}). "
+                "Returning best bound."
+            )
+            best = a if abs(fa) < abs(fb) else b
+            residual = min(abs(fa), abs(fb))
+            return EstimationResult(
+                value=best,
+                residual=residual,
+                iterations=iterations,
+                converged=residual <= self.tol,
+            )
+
+        for _ in range(self.maxiter):
+            mid = (a + b) / 2.0
+            fmid = objective(mid)
+            iterations += 1
+
+            if abs(fmid) <= self.tol or (b - a) / 2.0 <= self.tol:
+                return EstimationResult(
+                    value=mid,
+                    residual=abs(fmid),
+                    iterations=iterations,
+                    converged=True,
+                )
+
+            if fa * fmid < 0:
+                b = mid
+            else:
+                a = mid
+
+        mid = (a + b) / 2.0
+        return EstimationResult(
+            value=mid,
+            residual=abs(objective(mid)),
+            iterations=iterations + 1,
+            converged=False,
+        )
+
+
+class BrentMethod:
+    """Brent's method for root-finding (via scipy.optimize.brentq).
+
+    Requires that ``objective`` changes sign over ``bounds``.
+    Falls back to bisection if scipy is not available.
+
+    Args:
+        tol: Absolute tolerance on the root. Defaults to 1e-3.
+        maxiter: Maximum number of iterations. Defaults to 50.
+    """
+
+    def __init__(self, tol: float = 1e-3, maxiter: int = 50) -> None:
+        self.tol = tol
+        self.maxiter = maxiter
+
+    def solve(
+        self,
+        objective: Callable[[float], float],
+        bounds: tuple[float, float],
+    ) -> EstimationResult:
+        try:
+            from scipy.optimize import brentq  # type: ignore
+        except ImportError:
+            logger.warning(
+                "scipy not available, falling back to BisectionMethod."
+            )
+            fallback = BisectionMethod(tol=self.tol, maxiter=self.maxiter)
+            return fallback.solve(objective, bounds)
+
+        a, b = bounds
+
+        try:
+            root, result_info = brentq(
+                objective,
+                a,
+                b,
+                xtol=self.tol,
+                maxiter=self.maxiter,
+                full_output=True,
+            )
+            return EstimationResult(
+                value=root,
+                residual=abs(
+                    result_info.function_calls and abs(objective(root)) or 0.0
+                ),
+                iterations=result_info.iterations,
+                converged=result_info.converged,
+            )
+        except ValueError as e:
+            logger.warning(f"Brent's method failed: {e}. Trying bisection.")
+            fallback = BisectionMethod(tol=self.tol, maxiter=self.maxiter)
+            return fallback.solve(objective, bounds)
+
+
+class NewtonMethod:
+    """Newton-Raphson with finite-difference derivative.
+
+    Does not require bounds for the algorithm itself, but bounds are used
+    to clip iterates and provide the initial guess (midpoint).
+
+    Args:
+        tol: Absolute tolerance on the residual. Defaults to 1e-3.
+        maxiter: Maximum number of iterations. Defaults to 20.
+        dx: Step size for finite-difference derivative. Defaults to 1.0.
+    """
+
+    def __init__(
+        self,
+        tol: float = 1e-3,
+        maxiter: int = 20,
+        dx: float = 1.0,
+    ) -> None:
+        self.tol = tol
+        self.maxiter = maxiter
+        self.dx = dx
+
+    def solve(
+        self,
+        objective: Callable[[float], float],
+        bounds: tuple[float, float],
+    ) -> EstimationResult:
+        a, b = bounds
+        x = (a + b) / 2.0
+        iterations = 0
+
+        for _ in range(self.maxiter):
+            fx = objective(x)
+            iterations += 1
+
+            if abs(fx) <= self.tol:
+                return EstimationResult(
+                    value=x,
+                    residual=abs(fx),
+                    iterations=iterations,
+                    converged=True,
+                )
+
+            fx_dx = objective(x + self.dx)
+            iterations += 1
+            derivative = (fx_dx - fx) / self.dx
+
+            if abs(derivative) < 1e-12:
+                logger.warning(
+                    "Newton: near-zero derivative at x=%.6g, stopping.", x
+                )
+                return EstimationResult(
+                    value=x,
+                    residual=abs(fx),
+                    iterations=iterations,
+                    converged=False,
+                )
+
+            x_new = x - fx / derivative
+            # Clip to bounds
+            x_new = max(a, min(b, x_new))
+            x = x_new
+
+        fx = objective(x)
+        iterations += 1
+        return EstimationResult(
+            value=x,
+            residual=abs(fx),
+            iterations=iterations,
+            converged=abs(fx) <= self.tol,
+        )
diff --git a/src/mechaphlowers/core/models/estimation/result.py b/src/mechaphlowers/core/models/estimation/result.py
new file mode 100644
index 00000000..14580bcf
--- /dev/null
+++ b/src/mechaphlowers/core/models/estimation/result.py
@@ -0,0 +1,34 @@
+# Copyright (c) 2026, RTE (http://www.rte-france.com)
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+# SPDX-License-Identifier: MPL-2.0
+
+from __future__ import annotations
+
+from dataclasses import dataclass
+
+
+@dataclass(frozen=True)
+class EstimationResult:
+    """Result of an inverse estimation solve.
+
+    Attributes:
+        value: The estimated variable value that satisfies the target.
+        residual: Final residual ``|f(value) - target|``.
+        iterations: Number of algorithm iterations.
+        converged: Whether the algorithm converged within tolerance.
+    """
+
+    value: float
+    residual: float
+    iterations: int
+    converged: bool
+
+    def __repr__(self) -> str:
+        status = "converged" if self.converged else "NOT converged"
+        return (
+            f"EstimationResult(value={self.value:.6g}, "
+            f"residual={self.residual:.3e}, "
+            f"iterations={self.iterations}, {status})"
+        )
diff --git a/test/api/test_section_study.py b/test/api/test_section_study.py
index a4a905bd..fcd5aac5 100644
--- a/test/api/test_section_study.py
+++ b/test/api/test_section_study.py
@@ -166,45 +166,70 @@ def test_restore_notifies_position_engine(
 
 
 class TestSectionStudyLazyEngines:
-    def test_plot_engine_not_created_eagerly(
+    def test_engines_not_created_eagerly(
         self, balance_engine_base_test: BalanceEngine
-    ):
+    ) -> None:
         study = SectionStudy(
             cable_array=balance_engine_base_test.cable_array,
             section_array=balance_engine_base_test.section_array,
         )
         assert study._plot_engine is None
+        assert study._thermal_engine is None
+        assert study._guying is None
+        assert study._estimation_engine is None
 
     def test_plot_engine_created_on_access(
         self, balance_engine_base_test: BalanceEngine
-    ):
+    ) -> None:
         study = SectionStudy(
             cable_array=balance_engine_base_test.cable_array,
             section_array=balance_engine_base_test.section_array,
         )
-        _ = study.plot_engine
+        engine = study.plot_engine
         assert study._plot_engine is not None
 
+        # Second access returns same instance
+        assert study.plot_engine is engine
+
     def test_thermal_engine_created_on_access(
         self, balance_engine_base_test: BalanceEngine
-    ):
+    ) -> None:
         study = SectionStudy(
             cable_array=balance_engine_base_test.cable_array,
             section_array=balance_engine_base_test.section_array,
         )
-        _ = study.thermal_engine
+        engine = study.thermal_engine
         assert study._thermal_engine is not None
 
+        # Second access returns same instance
+        assert study.thermal_engine is engine
+
     def test_guying_created_on_access(
         self, balance_engine_base_test: BalanceEngine
-    ):
+    ) -> None:
         study = SectionStudy(
             cable_array=balance_engine_base_test.cable_array,
             section_array=balance_engine_base_test.section_array,
         )
-        _ = study.guying
+        guying = study.guying
         assert study._guying is not None
 
+        # Second access returns same instance
+        assert study.guying is guying
+
+    def test_estimation_engine_created_on_access(
+        self, balance_engine_base_test: BalanceEngine
+    ) -> None:
+        study = SectionStudy(
+            cable_array=balance_engine_base_test.cable_array,
+            section_array=balance_engine_base_test.section_array,
+        )
+        estimation_engine = study.estimation_engine
+        assert study._estimation_engine is not None
+
+        # Second access returns same instance
+        assert study.estimation_engine is estimation_engine
+
 
 class TestSectionStudySubEngines:
     def test_balance_engine_property(
diff --git a/test/core/models/test_estimation.py b/test/core/models/test_estimation.py
new file mode 100644
index 00000000..ecac0e5b
--- /dev/null
+++ b/test/core/models/test_estimation.py
@@ -0,0 +1,167 @@
+# Copyright (c) 2026, RTE (http://www.rte-france.com)
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+# SPDX-License-Identifier: MPL-2.0
+
+import numpy as np
+import pytest
+
+from mechaphlowers.api.section_study import SectionStudy
+from mechaphlowers.core.models.balance.engine import BalanceEngine
+from mechaphlowers.core.models.estimation import (
+    BisectionMethod,
+    BrentMethod,
+    EstimationEngine,
+    EstimationResult,
+    NewtonMethod,
+)
+
+
+@pytest.fixture
+def solved_study(balance_engine_base_test: BalanceEngine) -> SectionStudy:
+    """A SectionStudy that has already been solved via solve_adjustment."""
+    study = SectionStudy(
+        cable_array=balance_engine_base_test.cable_array,
+        section_array=balance_engine_base_test.section_array,
+    )
+    study.solve_adjustment()
+    return study
+
+
+class TestEstimationResult:
+    def test_repr_converged(self):
+        r = EstimationResult(
+            value=42.5, residual=1e-5, iterations=10, converged=True
+        )
+        assert "converged" in repr(r)
+        assert "42.5" in repr(r)
+
+    def test_repr_not_converged(self):
+        r = EstimationResult(
+            value=10.0, residual=5.0, iterations=50, converged=False
+        )
+        assert "NOT converged" in repr(r)
+
+
+class TestOptimizationMethods:
+    """Test the optimization methods on a simple analytical function."""
+
+    def _simple_objective(self, x: float) -> float:
+        """f(x) = x^2 - 4, root at x=2."""
+        return x**2 - 4
+
+    def test_bisection_converges(self):
+        method = BisectionMethod(tol=1e-6, maxiter=100)
+        result = method.solve(self._simple_objective, bounds=(0.0, 5.0))
+        assert result.converged
+        assert abs(result.value - 2.0) < 1e-5
+
+    def test_brent_converges(self):
+        method = BrentMethod(tol=1e-6, maxiter=100)
+        result = method.solve(self._simple_objective, bounds=(0.0, 5.0))
+        assert result.converged
+        assert abs(result.value - 2.0) < 1e-5
+
+    def test_newton_converges(self):
+        method = NewtonMethod(tol=1e-6, maxiter=50, dx=0.01)
+        result = method.solve(self._simple_objective, bounds=(0.0, 5.0))
+        assert result.converged
+        assert abs(result.value - 2.0) < 1e-4
+
+    def test_bisection_no_sign_change(self):
+        """When bounds don't bracket a root, bisection returns best bound."""
+        method = BisectionMethod(tol=1e-6, maxiter=50)
+        # f(3)=5, f(5)=21, both positive
+        result = method.solve(self._simple_objective, bounds=(3.0, 5.0))
+        assert not result.converged
+        assert result.value == 3.0
+
+
+class TestEstimationEngine:
+    def test_estimate_generic(self, solved_study: SectionStudy):
+        """Generic estimate() works with a simple function."""
+        engine = EstimationEngine(
+            solved_study, method=BisectionMethod(tol=1e-4)
+        )
+        result = engine.estimate(
+            objective=lambda x: x**2 - 9,
+            bounds=(0.0, 10.0),
+        )
+        assert result.converged
+        assert abs(result.value - 3.0) < 1e-3
+
+    def test_method_setter(self, solved_study: SectionStudy):
+        engine = EstimationEngine(solved_study)
+        assert isinstance(engine.method, BrentMethod)
+        engine.method = NewtonMethod(tol=1e-3)
+        assert isinstance(engine.method, NewtonMethod)
+
+    def test_estimate_temperature(self, solved_study: SectionStudy):
+        """Estimate temperature: solve at known temp, use that distance as target."""
+        # First, solve at a known temperature to get a reference distance
+        study = solved_study
+        temperature = 60.0
+        study.solve_change_state(new_temperature=temperature)
+        ref_distance = study.position_engine.point_distance(
+            span_index=0, point=np.array([250.0, 10.0, 20.0])
+        ).distance_3d
+
+        # Restore state
+        study.solve_adjustment()
+
+        # Now estimate what temperature gives that distance
+        engine = EstimationEngine(study, method=BrentMethod(tol=0.1))
+        result = engine.estimate_temperature(
+            span_index=0,
+            obstacle_point=np.array([250.0, 10.0, 20.0]),
+            target_distance=ref_distance,
+            bounds=(0.0, 120.0),
+        )
+        assert result.converged
+        assert abs(result.value - temperature) < 0.1
+
+    def test_estimate_wind(self, solved_study: SectionStudy):
+        """Estimate wind pressure using a known reference."""
+        study = solved_study
+        wind_pressure = 400.0
+        study.solve_change_state(wind_pressure=wind_pressure)
+        ref_distance = study.position_engine.point_distance(
+            span_index=0, point=np.array([250.0, 10.0, 20.0])
+        ).distance_3d
+
+        # Restore state
+        study.solve_adjustment()
+
+        engine = EstimationEngine(study, method=BrentMethod(tol=0.5))
+        result = engine.estimate_wind(
+            span_index=0,
+            obstacle_point=np.array([250.0, 10.0, 20.0]),
+            target_distance=ref_distance,
+            bounds=(0.0, 1000.0),
+        )
+        assert result.converged
+        assert abs(result.value - wind_pressure) < 0.1
+
+    def test_state_preserved_after_estimation(
+        self, solved_study: SectionStudy
+    ):
+        """Engine state is unchanged after estimation."""
+        study = solved_study
+        memento_before = study.save_state()
+
+        engine = EstimationEngine(
+            study, method=BisectionMethod(tol=1.0, maxiter=5)
+        )
+        engine.estimate_temperature(
+            span_index=0,
+            obstacle_point=np.array([250.0, 10.0, 20.0]),
+            target_distance=50.0,
+            bounds=(0.0, 100.0),
+        )
+
+        # State should be unchanged
+        np.testing.assert_array_almost_equal(
+            study.balance_engine.balance_model.nodes.dxdydz,
+            memento_before.dxdydz,
+        )