diff --git a/.claude/agents/reviewer.md b/.claude/agents/reviewer.md
index 8179420e..fbabf5c1 100644
--- a/.claude/agents/reviewer.md
+++ b/.claude/agents/reviewer.md
@@ -48,7 +48,7 @@ End with totals + verdict: APPROVE / REQUEST CHANGES.
 - [ ] SOLID: one concern per class, constructor injection, depend on abstractions, no duplicated logic.
 
 **Interfaces & errors (WARNING)**
-- [ ] `virtual ~I() = default` — never `= 0`. No exceptions — `std::optional`/status enums.
+- [ ] Interfaces: `protected:` non-virtual `~I() = default` — no `virtual` or `= 0` destructor, no ownership or deletion through an interface pointer. No exceptions — `std::optional`/status enums.
 
 **Testing (WARNING)**
 - [ ] `TEST_F` on `float` — no `TYPED_TEST`, no multi-type, never plain `TEST()`.
diff --git a/.github/agents/reviewer.agent.md b/.github/agents/reviewer.agent.md
index 5b7ac011..a7cbe005 100644
--- a/.github/agents/reviewer.agent.md
+++ b/.github/agents/reviewer.agent.md
@@ -54,7 +54,7 @@ End with totals + verdict: APPROVE / REQUEST CHANGES.
 - [ ] SOLID: one concern per class, constructor injection, depend on abstractions, no duplicated logic.
 
 **Interfaces & errors (WARNING)**
-- [ ] `virtual ~I() = default` — never `= 0`. No exceptions — `std::optional`/status enums.
+- [ ] Interfaces: `protected:` non-virtual `~I() = default` — no `virtual` or `= 0` destructor, no ownership or deletion through an interface pointer. No exceptions — `std::optional`/status enums.
 
 **Testing (WARNING)**
 - [ ] `TEST_F` on `float` — no `TYPED_TEST`, no multi-type, never plain `TEST()`.
diff --git a/.github/copilot-instructions.md b/.github/copilot-instructions.md
index f2a71873..dea28ed9 100644
--- a/.github/copilot-instructions.md
+++ b/.github/copilot-instructions.md
@@ -9,5 +9,5 @@ Essentials (full detail in AGENTS.md):
 - **Embedded** — no `#pragma GCC optimize` / `optimize` attribute (they block inlining); the consumer sets optimization flags per TU; `OPTIMIZE_FOR_SPEED` (forced inlining) on hot paths.
 - **No comments** (except license/NOLINT). Allman braces, brace-init, PascalCase/camelCase.
 - **Tests** — `TEST_F` on `float`, `StrictMock` only, never plain `TEST()`, no redundant cases.
-- **No exceptions** — `std::optional`/status enums; interfaces `virtual ~I() = default`.
+- **No exceptions** — `std::optional`/status enums; interfaces get a `protected:` non-virtual `~I() = default`, never a virtual destructor.
 - **Be terse** — minimal prose; report file paths + pass/fail.
diff --git a/.github/instructions/numerical-cpp.instructions.md b/.github/instructions/numerical-cpp.instructions.md
index 7bd08028..7ac6ba25 100644
--- a/.github/instructions/numerical-cpp.instructions.md
+++ b/.github/instructions/numerical-cpp.instructions.md
@@ -68,7 +68,7 @@ Use `math::IsFinite` for finiteness checks: a consumer may still build with `-ff
 - SOLID principles — constructor injection, depend on abstractions
 - DRY — extract shared logic into helpers or templates
 - RAII for resource management
-- **No pure virtual destructors**: Never declare `virtual ~Interface() = 0`. Pure virtual destructors force the compiler to emit a thunk and a separate body, wasting memory on embedded targets. Use `virtual ~Interface() = default` if a virtual destructor is required, or omit it entirely when the class is never deleted through a base pointer.
+- **Protected, non-virtual interface destructors**: declare `~Interface() = default` in a `protected:` section, never `virtual ~Interface()` (defaulted or `= 0`). A virtual destructor puts a deleting destructor in every implementation's vtable, which links `operator delete` (and with it `malloc`/`free`) into heap-free firmware, and under multiple inheritance adds a thunk per base. The protected destructor makes deletion through an interface pointer a compile error. Concrete classes declare no destructor unless they need one, and callers hold the concrete type rather than a `std::unique_ptr` to the interface.
 - No virtual calls in ISR-callable or real-time critical paths
 
 ## Documentation — MANDATORY
diff --git a/AGENTS.md b/AGENTS.md
index 08f5c37e..6db01b3f 100644
--- a/AGENTS.md
+++ b/AGENTS.md
@@ -44,7 +44,8 @@ resource-constrained embedded systems. Real-time, deterministic, no heap.
 
 ## Interfaces & errors
 
-- Interfaces = pure virtual classes; `virtual ~I() = default` — **never** `= 0` destructors.
+- Interfaces = pure virtual classes with a `protected:`, non-virtual `~I() = default` — **never** a `virtual` or `= 0` destructor. A virtual one puts a deleting destructor in every vtable, which links `operator delete`/`free` into heap-free firmware.
+- Concrete classes declare no destructor unless they need one. Never own or delete through an interface pointer (no `std::unique_ptr<I>`): hold the concrete type, pass `I&`.
 - No exceptions. Use `std::optional<T>` or status enums. `assert`/`really_assert` for preconditions.
 - SOLID + DIP: constructor injection; depend on abstractions. RAII. No virtual calls in real-time paths.
 
diff --git a/CLAUDE.md b/CLAUDE.md
index 994b533b..1ff8398d 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -9,5 +9,5 @@ Essentials (full detail in AGENTS.md):
 - **Embedded** — no `#pragma GCC optimize` / `optimize` attribute (they block inlining); the consumer sets optimization flags per TU; `OPTIMIZE_FOR_SPEED` (forced inlining) on hot paths.
 - **No comments** (except license/NOLINT). Allman braces, brace-init, PascalCase/camelCase.
 - **Tests** — `TEST_F` on `float`, `StrictMock` only, never plain `TEST()`, no redundant cases.
-- **No exceptions** — `std::optional`/status enums; interfaces `virtual ~I() = default`.
+- **No exceptions** — `std::optional`/status enums; interfaces get a `protected:` non-virtual `~I() = default`, never a virtual destructor.
 - **Be terse** — minimal prose; report file paths + pass/fail.
diff --git a/numerical/analysis/FastFourierTransform.hpp b/numerical/analysis/FastFourierTransform.hpp
index 6d163c26..4e9bef5e 100644
--- a/numerical/analysis/FastFourierTransform.hpp
+++ b/numerical/analysis/FastFourierTransform.hpp
@@ -9,8 +9,10 @@ namespace analysis
     class TwiddleFactors
     {
     public:
-        virtual ~TwiddleFactors() = default;
         virtual math::Complex<QNumberType>& operator[](std::size_t n) = 0;
+
+    protected:
+        ~TwiddleFactors() = default;
     };
 
     template<typename QNumberType>
@@ -21,8 +23,6 @@ namespace analysis
             "FastFourierTransform can only be instantiated with math::QNumber types.");
 
     public:
-        virtual ~FastFourierTransform() = default;
-
         static constexpr std::size_t Log2(std::size_t n)
         {
             std::size_t result = 0;
@@ -41,6 +41,8 @@ namespace analysis
         virtual VectorReal& Inverse(VectorComplex& input) = 0;
 
     protected:
+        ~FastFourierTransform() = default;
+
         template<typename T>
         inline T BitReverse(T x, T logr_n, T bits, T radix)
         {
diff --git a/numerical/analysis/windowing/Windowing.hpp b/numerical/analysis/windowing/Windowing.hpp
index ef825db9..153834f1 100644
--- a/numerical/analysis/windowing/Windowing.hpp
+++ b/numerical/analysis/windowing/Windowing.hpp
@@ -14,9 +14,11 @@ namespace windowing
             "Window can only be instantiated with math::QNumber types.");
 
     public:
-        virtual ~Window() = default;
         virtual QNumberType operator()(std::size_t n, std::size_t order) = 0;
         virtual QNumberType Power(std::size_t order) = 0;
+
+    protected:
+        ~Window() = default;
     };
 
     namespace detail
diff --git a/numerical/analysis/windowing/test/TestWindowing.cpp b/numerical/analysis/windowing/test/TestWindowing.cpp
index 68197bc2..8ee79ff3 100644
--- a/numerical/analysis/windowing/test/TestWindowing.cpp
+++ b/numerical/analysis/windowing/test/TestWindowing.cpp
@@ -11,7 +11,6 @@ namespace
     {
     protected:
         static constexpr float kEpsilon = 1e-3f;
-        std::unique_ptr<windowing::Window<QNumType>> window;
     };
 
     using WindowingTypes = ::testing::Types<math::Q31, math::Q15, float>;
@@ -20,7 +19,7 @@ namespace
 
 TYPED_TEST(WindowingTest, HammingWindowValues)
 {
-    this->window = std::make_unique<windowing::HammingWindow<TypeParam>>();
+    windowing::HammingWindow<TypeParam> window;
 
     struct TestCase
     {
@@ -36,7 +35,7 @@ TYPED_TEST(WindowingTest, HammingWindowValues)
 
     for (const auto& testCase : testCases)
     {
-        TypeParam result = (*this->window)(testCase.n, testCase.order);
+        TypeParam result = window(testCase.n, testCase.order);
         EXPECT_NEAR(math::ToFloat(result), testCase.expected, this->kEpsilon)
             << "Failed for n=" << testCase.n
             << ", order=" << testCase.order;
@@ -45,7 +44,7 @@ TYPED_TEST(WindowingTest, HammingWindowValues)
 
 TYPED_TEST(WindowingTest, HanningWindowValues)
 {
-    this->window = std::make_unique<windowing::HanningWindow<TypeParam>>();
+    windowing::HanningWindow<TypeParam> window;
 
     struct TestCase
     {
@@ -61,7 +60,7 @@ TYPED_TEST(WindowingTest, HanningWindowValues)
 
     for (const auto& testCase : testCases)
     {
-        TypeParam result = (*this->window)(testCase.n, testCase.order);
+        TypeParam result = window(testCase.n, testCase.order);
         EXPECT_NEAR(math::ToFloat(result), testCase.expected, this->kEpsilon)
             << "Failed for n=" << testCase.n
             << ", order=" << testCase.order;
@@ -70,7 +69,7 @@ TYPED_TEST(WindowingTest, HanningWindowValues)
 
 TYPED_TEST(WindowingTest, BlackmanWindowValues)
 {
-    this->window = std::make_unique<windowing::BlackmanWindow<TypeParam>>();
+    windowing::BlackmanWindow<TypeParam> window;
 
     struct TestCase
     {
@@ -85,7 +84,7 @@ TYPED_TEST(WindowingTest, BlackmanWindowValues)
 
     for (const auto& testCase : testCases)
     {
-        TypeParam result = (*this->window)(testCase.n, testCase.order);
+        TypeParam result = window(testCase.n, testCase.order);
         EXPECT_NEAR(math::ToFloat(result), testCase.expected, this->kEpsilon)
             << "Failed for n=" << testCase.n
             << ", order=" << testCase.order;
@@ -94,11 +93,11 @@ TYPED_TEST(WindowingTest, BlackmanWindowValues)
 
 TYPED_TEST(WindowingTest, RectangularWindowValues)
 {
-    this->window = std::make_unique<windowing::RectangularWindow<TypeParam>>();
+    windowing::RectangularWindow<TypeParam> window;
 
     for (std::size_t n = 0; n <= 8; n += 2)
     {
-        TypeParam result = (*this->window)(n, 8);
+        TypeParam result = window(n, 8);
         EXPECT_NEAR(math::ToFloat(result), 0.999f, this->kEpsilon)
             << "Failed for n=" << n;
     }
@@ -106,14 +105,14 @@ TYPED_TEST(WindowingTest, RectangularWindowValues)
 
 TYPED_TEST(WindowingTest, WindowSymmetry)
 {
-    std::vector<std::unique_ptr<windowing::Window<TypeParam>>> windows;
-    windows.push_back(std::make_unique<windowing::HammingWindow<TypeParam>>());
-    windows.push_back(std::make_unique<windowing::HanningWindow<TypeParam>>());
-    windows.push_back(std::make_unique<windowing::BlackmanWindow<TypeParam>>());
+    windowing::HammingWindow<TypeParam> hamming;
+    windowing::HanningWindow<TypeParam> hanning;
+    windowing::BlackmanWindow<TypeParam> blackman;
+    std::array<windowing::Window<TypeParam>*, 3> windows{ &hamming, &hanning, &blackman };
 
     const std::size_t order = 16;
 
-    for (const auto& window : windows)
+    for (auto* window : windows)
     {
         for (std::size_t n = 0; n < order / 2; ++n)
         {
diff --git a/numerical/controllers/implementations/Feedforward2Dof.hpp b/numerical/controllers/implementations/Feedforward2Dof.hpp
index 365001b0..5da2f209 100644
--- a/numerical/controllers/implementations/Feedforward2Dof.hpp
+++ b/numerical/controllers/implementations/Feedforward2Dof.hpp
@@ -12,9 +12,11 @@ namespace controllers
     public:
         static_assert(std::is_floating_point_v<T>, "Feedforward supports floating-point types");
 
-        virtual ~Feedforward() = default;
         virtual T Evaluate(T reference) = 0;
         virtual void Reset() = 0;
+
+    protected:
+        ~Feedforward() = default;
     };
 
     template<typename T>
@@ -23,9 +25,11 @@ namespace controllers
     public:
         static_assert(std::is_floating_point_v<T>, "FeedbackLaw supports floating-point types");
 
-        virtual ~FeedbackLaw() = default;
         virtual T Process(T error) = 0;
         virtual void Reset() = 0;
+
+    protected:
+        ~FeedbackLaw() = default;
     };
 
     template<typename T>
diff --git a/numerical/controllers/interfaces/OutputFeedbackController.hpp b/numerical/controllers/interfaces/OutputFeedbackController.hpp
index f911b472..c7c3516b 100644
--- a/numerical/controllers/interfaces/OutputFeedbackController.hpp
+++ b/numerical/controllers/interfaces/OutputFeedbackController.hpp
@@ -21,7 +21,9 @@ namespace controllers
         using MeasurementVector = math::Vector<T, MeasurementSize>;
         using InputVector = math::Vector<T, InputSize>;
 
-        virtual ~OutputFeedbackController() = default;
         virtual InputVector ComputeControl(const MeasurementVector& measurement) = 0;
+
+    protected:
+        ~OutputFeedbackController() = default;
     };
 }
diff --git a/numerical/controllers/interfaces/PidController.hpp b/numerical/controllers/interfaces/PidController.hpp
index b0cb7656..656d031b 100644
--- a/numerical/controllers/interfaces/PidController.hpp
+++ b/numerical/controllers/interfaces/PidController.hpp
@@ -26,10 +26,12 @@ namespace controllers
             "Pid can only be instantiated with math::QNumber types.");
 
     public:
-        virtual ~AsynchronousPidController() = default;
         virtual void Reset() = 0;
         virtual void SetTunings(PidTunings<QNumberType> tunings) = 0;
         virtual void SetPoint(QNumberType setPoint) = 0;
+
+    protected:
+        ~AsynchronousPidController() = default;
     };
 
     template<typename QNumberType>
@@ -38,5 +40,8 @@ namespace controllers
     {
     public:
         virtual QNumberType Process(QNumberType processVariable) = 0;
+
+    protected:
+        ~SynchronousPidController() = default;
     };
 }
diff --git a/numerical/controllers/interfaces/PidDriver.hpp b/numerical/controllers/interfaces/PidDriver.hpp
index f58bd11c..d58bbff2 100644
--- a/numerical/controllers/interfaces/PidDriver.hpp
+++ b/numerical/controllers/interfaces/PidDriver.hpp
@@ -9,10 +9,12 @@ namespace controllers
     class PidDriver
     {
     public:
-        virtual ~PidDriver() = default;
         virtual void Read(const infra::Function<void(QNumberType)>& onDone) = 0;
         virtual void ControlAction(QNumberType) = 0;
         virtual void Start(std::chrono::system_clock::duration sampleTime) = 0;
         virtual void Stop() = 0;
+
+    protected:
+        ~PidDriver() = default;
     };
 }
diff --git a/numerical/controllers/interfaces/StateFeedbackController.hpp b/numerical/controllers/interfaces/StateFeedbackController.hpp
index f8f8dab3..74f3eea4 100644
--- a/numerical/controllers/interfaces/StateFeedbackController.hpp
+++ b/numerical/controllers/interfaces/StateFeedbackController.hpp
@@ -19,7 +19,9 @@ namespace controllers
         using StateVector = math::Vector<T, StateSize>;
         using InputVector = math::Vector<T, InputSize>;
 
-        virtual ~StateFeedbackController() = default;
         virtual InputVector ComputeControl(const StateVector& state) = 0;
+
+    protected:
+        ~StateFeedbackController() = default;
     };
 }
diff --git a/numerical/estimators/Estimator.hpp b/numerical/estimators/Estimator.hpp
index f3c669fb..caebd7c5 100644
--- a/numerical/estimators/Estimator.hpp
+++ b/numerical/estimators/Estimator.hpp
@@ -11,8 +11,6 @@ namespace estimators
             "OfflineEstimator only supports float or QNumber types");
 
     public:
-        virtual ~OfflineEstimator() = default;
-
         using CoefficientsMatrix = math::Matrix<T, Features + 1, 1>;
         using DesignMatrix = math::Matrix<T, Features + 1, Features + 1>;
         using InputMatrix = math::Matrix<T, Features, 1>;
@@ -20,6 +18,9 @@ namespace estimators
         virtual bool Fit(const math::Matrix<T, Samples, Features>& X, const math::Matrix<T, Samples, 1>& y) = 0;
         virtual T Predict(const InputMatrix& X) const = 0;
         virtual const CoefficientsMatrix& Coefficients() const = 0;
+
+    protected:
+        ~OfflineEstimator() = default;
     };
 
     template<typename T, std::size_t Features>
@@ -29,8 +30,6 @@ namespace estimators
             "OnlineEstimator only supports float or QNumber types");
 
     public:
-        virtual ~OnlineEstimator() = default;
-
         using CoefficientsMatrix = math::Matrix<T, Features, 1>;
         using DesignMatrix = math::Matrix<T, Features, Features>;
         using InputMatrix = math::Matrix<T, Features, 1>;
@@ -44,5 +43,8 @@ namespace estimators
 
         virtual EstimationMetrics Update(const InputMatrix& x, const math::Matrix<T, 1, 1>& y) = 0;
         virtual const CoefficientsMatrix& Coefficients() const = 0;
+
+    protected:
+        ~OnlineEstimator() = default;
     };
 }
diff --git a/numerical/nonlinear_control/BacksteppingControl.hpp b/numerical/nonlinear_control/BacksteppingControl.hpp
index 3c33f5b9..c831ee57 100644
--- a/numerical/nonlinear_control/BacksteppingControl.hpp
+++ b/numerical/nonlinear_control/BacksteppingControl.hpp
@@ -18,11 +18,12 @@ namespace nonlinear_control
     public:
         using StateVector = std::array<T, Order>;
 
-        virtual ~StrictFeedbackModel() = default;
-
         virtual T Drift(std::size_t i, const StateVector& x) const = 0;
         virtual T Gain(std::size_t i, const StateVector& x) const = 0;
         virtual T VirtualDerivative(std::size_t i, const StateVector& x, T alpha) const = 0;
+
+    protected:
+        ~StrictFeedbackModel() = default;
     };
 
     template<typename T, std::size_t Order>
diff --git a/numerical/nonlinear_control/FeedbackLinearization.hpp b/numerical/nonlinear_control/FeedbackLinearization.hpp
index 44785605..3b6da8a8 100644
--- a/numerical/nonlinear_control/FeedbackLinearization.hpp
+++ b/numerical/nonlinear_control/FeedbackLinearization.hpp
@@ -18,10 +18,11 @@ namespace nonlinear_control
         using StateVector = math::Vector<T, Dim>;
         using DecouplingMatrix = math::SquareMatrix<T, Dim>;
 
-        virtual ~ControlAffineModel() = default;
-
         [[nodiscard]] virtual DecouplingMatrix DecouplingMatrixAt(const StateVector& x) const = 0;
         [[nodiscard]] virtual StateVector DriftTerm(const StateVector& x) const = 0;
+
+    protected:
+        ~ControlAffineModel() = default;
     };
 
     template<typename T, std::size_t Dim>
diff --git a/numerical/optimization/BayesianOptimization.hpp b/numerical/optimization/BayesianOptimization.hpp
index b54aa85e..b0d2140b 100644
--- a/numerical/optimization/BayesianOptimization.hpp
+++ b/numerical/optimization/BayesianOptimization.hpp
@@ -52,8 +52,10 @@ namespace optimization
         class BlackBoxObjective
         {
         public:
-            virtual ~BlackBoxObjective() = default;
             virtual float Evaluate(const ParameterVector& params) = 0;
+
+        protected:
+            ~BlackBoxObjective() = default;
         };
 
         BayesianOptimization(const BoundsArray& bounds,
diff --git a/numerical/optimization/DampingPolicy.hpp b/numerical/optimization/DampingPolicy.hpp
index 9fa9a597..2e9459e5 100644
--- a/numerical/optimization/DampingPolicy.hpp
+++ b/numerical/optimization/DampingPolicy.hpp
@@ -12,9 +12,10 @@ namespace optimization
         static_assert(std::is_floating_point_v<T>, "DampingPolicy supports floating-point types only");
 
     public:
-        virtual ~DampingPolicy() = default;
-
         [[nodiscard]] virtual T DampingSquared(const math::Matrix<T, Residuals, Parameters>& jacobian) const = 0;
+
+    protected:
+        ~DampingPolicy() = default;
     };
 
     template<typename T, std::size_t Residuals, std::size_t Parameters>
diff --git a/numerical/optimization/LeastSquaresProblem.hpp b/numerical/optimization/LeastSquaresProblem.hpp
index 2afb9fdc..ed6cad0c 100644
--- a/numerical/optimization/LeastSquaresProblem.hpp
+++ b/numerical/optimization/LeastSquaresProblem.hpp
@@ -12,10 +12,11 @@ namespace optimization
         static_assert(std::is_floating_point_v<T>, "LeastSquaresProblem supports floating-point types only");
 
     public:
-        virtual ~LeastSquaresProblem() = default;
-
         [[nodiscard]] virtual math::Vector<T, Residuals> Residual(const State& x) const = 0;
         [[nodiscard]] virtual math::Matrix<T, Residuals, Parameters> Jacobian(const State& x) const = 0;
         [[nodiscard]] virtual State Retract(const State& x, const math::Vector<T, Parameters>& delta) const = 0;
+
+    protected:
+        ~LeastSquaresProblem() = default;
     };
 }
diff --git a/numerical/optimization/ObjectiveFunction.hpp b/numerical/optimization/ObjectiveFunction.hpp
index 2d3152f6..0dbe478c 100644
--- a/numerical/optimization/ObjectiveFunction.hpp
+++ b/numerical/optimization/ObjectiveFunction.hpp
@@ -15,6 +15,8 @@ namespace optimization
 
         virtual QNumberType Cost(const Vector& parameters) = 0;
         virtual Vector Gradient(const Vector& parameters) = 0;
-        virtual ~ObjectiveFunction() = default;
+
+    protected:
+        ~ObjectiveFunction() = default;
     };
 }
diff --git a/numerical/optimization/Optimizer.hpp b/numerical/optimization/Optimizer.hpp
index 3e87d09d..8def38ae 100644
--- a/numerical/optimization/Optimizer.hpp
+++ b/numerical/optimization/Optimizer.hpp
@@ -12,8 +12,6 @@ namespace optimization
             "Optimizer can only be instantiated with math::QNumber types.");
 
     public:
-        virtual ~Optimizer() = default;
-
         using Vector = math::Vector<QNumberType, NumberOfFeatures>;
 
         struct Result
@@ -30,5 +28,8 @@ namespace optimization
         };
 
         virtual const Result& Minimize(const Vector& initialGuess, ObjectiveFunction<QNumberType, NumberOfFeatures>& objective) = 0;
+
+    protected:
+        ~Optimizer() = default;
     };
 }
diff --git a/numerical/optimization/StepOptimizer.hpp b/numerical/optimization/StepOptimizer.hpp
index 9859b4c7..0ca871a1 100644
--- a/numerical/optimization/StepOptimizer.hpp
+++ b/numerical/optimization/StepOptimizer.hpp
@@ -10,11 +10,12 @@ namespace optimization
         static_assert(std::is_floating_point_v<T>, "StepOptimizer supports floating-point types only");
 
     public:
-        virtual ~StepOptimizer() = default;
-
         using Vector = math::Vector<T, N>;
 
         virtual void Step(Vector& theta, const Vector& gradient) = 0;
         virtual void Reset() = 0;
+
+    protected:
+        ~StepOptimizer() = default;
     };
 }
diff --git a/numerical/regularization/Regularization.hpp b/numerical/regularization/Regularization.hpp
index 8a41fba5..5bee29a4 100644
--- a/numerical/regularization/Regularization.hpp
+++ b/numerical/regularization/Regularization.hpp
@@ -14,6 +14,8 @@ namespace regularization
         using Vector = math::Vector<T, Size>;
         virtual T Calculate(const Vector& parameters) const = 0;
         virtual Vector Gradient(const Vector& parameters) const = 0;
-        virtual ~Regularization() = default;
+
+    protected:
+        ~Regularization() = default;
     };
 }
diff --git a/numerical/solvers/OdeSystem.hpp b/numerical/solvers/OdeSystem.hpp
index 1a21d3ff..0d386aa0 100644
--- a/numerical/solvers/OdeSystem.hpp
+++ b/numerical/solvers/OdeSystem.hpp
@@ -12,11 +12,13 @@ namespace solvers
         static_assert(std::is_floating_point_v<T>, "OdeSystem supports floating-point types only");
 
     public:
-        virtual ~OdeSystem() = default;
         [[nodiscard]] virtual math::Vector<T, StateSize> Derivative(
             const math::Vector<T, StateSize>& x,
             const math::Vector<T, InputSize>& u,
             T t) = 0;
+
+    protected:
+        ~OdeSystem() = default;
     };
 
     template<typename T, std::size_t StateSize>
@@ -25,9 +27,11 @@ namespace solvers
         static_assert(std::is_floating_point_v<T>, "OdeSystem supports floating-point types only");
 
     public:
-        virtual ~OdeSystem() = default;
         [[nodiscard]] virtual math::Vector<T, StateSize> Derivative(
             const math::Vector<T, StateSize>& x,
             T t) = 0;
+
+    protected:
+        ~OdeSystem() = default;
     };
 }
diff --git a/numerical/solvers/SecondOrderSystem.hpp b/numerical/solvers/SecondOrderSystem.hpp
index 070d85f6..de0bcbd5 100644
--- a/numerical/solvers/SecondOrderSystem.hpp
+++ b/numerical/solvers/SecondOrderSystem.hpp
@@ -12,10 +12,12 @@ namespace solvers
         static_assert(std::is_floating_point_v<T>, "SecondOrderSystem supports floating-point types only");
 
     public:
-        virtual ~SecondOrderSystem() = default;
         [[nodiscard]] virtual math::Vector<T, Dimension> Acceleration(
             const math::Vector<T, Dimension>& q,
             const math::Vector<T, Dimension>& qDot,
             const math::Vector<T, InputSize>& u) = 0;
+
+    protected:
+        ~SecondOrderSystem() = default;
     };
 }
diff --git a/numerical/solvers/Solver.hpp b/numerical/solvers/Solver.hpp
index 41c7c18e..a730ce12 100644
--- a/numerical/solvers/Solver.hpp
+++ b/numerical/solvers/Solver.hpp
@@ -13,12 +13,13 @@ namespace solvers
             "Solver dimensions must be positive");
 
     public:
-        virtual ~Solver() = default;
-
         using SolutionVector = math::Vector<T, N>;
         using InputMatrix = math::Matrix<T, N, N>;
         using InputVector = math::Vector<T, N>;
 
         virtual SolutionVector Solve(const InputMatrix& a, const InputVector& b) = 0;
+
+    protected:
+        ~Solver() = default;
     };
 }
diff --git a/simulator/analysis/FastFourierTransform/application/FftSimulator.cpp b/simulator/analysis/FastFourierTransform/application/FftSimulator.cpp
index ceeb313e..3a97b0bb 100644
--- a/simulator/analysis/FastFourierTransform/application/FftSimulator.cpp
+++ b/simulator/analysis/FastFourierTransform/application/FftSimulator.cpp
@@ -1,40 +1,43 @@
 #include "simulator/analysis/FastFourierTransform/application/FftSimulator.hpp"
 #include "numerical/analysis/FastFourierTransformRadix2Impl.hpp"
-#include "numerical/analysis/windowing/Windowing.hpp"
 #include "numerical/analysis/TwiddleFactorsTable.hpp"
+#include "numerical/analysis/windowing/Windowing.hpp"
 #include <algorithm>
 #include <cmath>
-#include <memory>
 #include <stdexcept>
 
 namespace simulator::analysis
 {
     namespace
     {
-        std::vector<float> ApplyWindow(const std::vector<float>& signal, WindowType windowType)
+        windowing::Window<float>& SelectWindow(WindowType windowType)
         {
-            std::unique_ptr<windowing::Window<float>> window;
+            static windowing::HammingWindow<float> hamming;
+            static windowing::HanningWindow<float> hanning;
+            static windowing::BlackmanWindow<float> blackman;
+            static windowing::RectangularWindow<float> rectangular;
 
             switch (windowType)
             {
                 case WindowType::Hamming:
-                    window = std::make_unique<windowing::HammingWindow<float>>();
-                    break;
+                    return hamming;
                 case WindowType::Hanning:
-                    window = std::make_unique<windowing::HanningWindow<float>>();
-                    break;
+                    return hanning;
                 case WindowType::Blackman:
-                    window = std::make_unique<windowing::BlackmanWindow<float>>();
-                    break;
+                    return blackman;
                 case WindowType::Rectangular:
                 default:
-                    window = std::make_unique<windowing::RectangularWindow<float>>();
-                    break;
+                    return rectangular;
             }
+        }
+
+        std::vector<float> ApplyWindow(const std::vector<float>& signal, WindowType windowType)
+        {
+            windowing::Window<float>& window{ SelectWindow(windowType) };
 
             std::vector<float> windowed(signal.size());
             for (std::size_t i = 0; i < signal.size(); ++i)
-                windowed[i] = signal[i] * (*window)(i, signal.size());
+                windowed[i] = signal[i] * window(i, signal.size());
 
             return windowed;
         }
diff --git a/simulator/analysis/PowerDensitySpectrum/application/PsdSimulator.cpp b/simulator/analysis/PowerDensitySpectrum/application/PsdSimulator.cpp
index 6319b764..28e78e61 100644
--- a/simulator/analysis/PowerDensitySpectrum/application/PsdSimulator.cpp
+++ b/simulator/analysis/PowerDensitySpectrum/application/PsdSimulator.cpp
@@ -1,17 +1,37 @@
 #include "simulator/analysis/PowerDensitySpectrum/application/PsdSimulator.hpp"
 #include "numerical/analysis/FastFourierTransformRadix2Impl.hpp"
 #include "numerical/analysis/PowerDensitySpectrum.hpp"
-#include "numerical/analysis/windowing/Windowing.hpp"
 #include "numerical/analysis/TwiddleFactorsTable.hpp"
+#include "numerical/analysis/windowing/Windowing.hpp"
 #include <algorithm>
 #include <cmath>
-#include <memory>
 #include <stdexcept>
 
 namespace simulator::analysis::psd
 {
     namespace
     {
+        windowing::Window<float>& SelectWindow(WindowType windowType)
+        {
+            static windowing::HammingWindow<float> hamming;
+            static windowing::HanningWindow<float> hanning;
+            static windowing::BlackmanWindow<float> blackman;
+            static windowing::RectangularWindow<float> rectangular;
+
+            switch (windowType)
+            {
+                case WindowType::Hamming:
+                    return hamming;
+                case WindowType::Hanning:
+                    return hanning;
+                case WindowType::Blackman:
+                    return blackman;
+                case WindowType::Rectangular:
+                default:
+                    return rectangular;
+            }
+        }
+
         template<std::size_t SegmentSize, std::size_t Overlap>
         PsdResult ComputeForSizeAndOverlap(float sampleRateHz, std::size_t inputSize,
             const std::vector<float>& signal, windowing::Window<float>& window)
@@ -94,36 +114,20 @@ namespace simulator::analysis::psd
         auto signal = utils::SignalGenerator::Generate(configuration.inputSize, configuration.sampleRateHz,
             configuration.signalComponents, configuration.noise);
 
-        std::unique_ptr<windowing::Window<float>> window;
-        switch (configuration.windowType)
-        {
-            case WindowType::Hamming:
-                window = std::make_unique<windowing::HammingWindow<float>>();
-                break;
-            case WindowType::Hanning:
-                window = std::make_unique<windowing::HanningWindow<float>>();
-                break;
-            case WindowType::Blackman:
-                window = std::make_unique<windowing::BlackmanWindow<float>>();
-                break;
-            case WindowType::Rectangular:
-            default:
-                window = std::make_unique<windowing::RectangularWindow<float>>();
-                break;
-        }
+        windowing::Window<float>& window{ SelectWindow(configuration.windowType) };
 
         switch (configuration.segmentSize)
         {
             case 64:
-                return DispatchOverlap<64>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, *window);
+                return DispatchOverlap<64>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, window);
             case 128:
-                return DispatchOverlap<128>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, *window);
+                return DispatchOverlap<128>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, window);
             case 256:
-                return DispatchOverlap<256>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, *window);
+                return DispatchOverlap<256>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, window);
             case 512:
-                return DispatchOverlap<512>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, *window);
+                return DispatchOverlap<512>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, window);
             case 1024:
-                return DispatchOverlap<1024>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, *window);
+                return DispatchOverlap<1024>(configuration.overlapPercent, configuration.sampleRateHz, configuration.inputSize, signal, window);
             default:
                 throw std::invalid_argument("Unsupported segment size. Must be a power of 2 between 64 and 1024.");
         }