feat(cuk): make maximumDutyCycle configurable, replace hardcoded 0.95

Cuk::calculate_duty_cycle previously hard-coded a 0.95 ceiling and
threw on D >= 0.95. That was already loud (no silent clamp) but
inflexible -- callers with stricter controller constraints (e.g.
maxD = 0.7) had no way to express that and got false-positive
"working" designs that violated their actual control range.

Adds std::optional<double> maximumDutyCycle = 0.95 field with
get/set accessors, and threads it through calculate_duty_cycle as a
new optional 6th parameter (default 0.95 preserves backward
compatibility for the existing tests in TestCuk.cpp). All 9 internal
call sites updated to pass maximumDutyCycle.value_or(0.95). The
inequality also gains the standard 1 % rounding tolerance to match
the rest of the converter family.

Tests: all 48 [cuk-topology] cases pass (324 assertions).

Phase C/5 of 9. Same pattern as Flyback (04272d7), forward family
(683e731), Buck (2c9300c), Boost (96fdb52), IsolatedBuck
(703bc80), IsolatedBuckBoost (a158d54).

Author: AlfVII

src/converter_models/Cuk.cpp

@@ -16,7 +16,7 @@ namespace OpenMagnetics {
     //   Static analytical helpers (CUK_PLAN.md §2.13, §4.1–§4.4)
     // ============================================================
 
-    double Cuk::calculate_duty_cycle(double inputVoltage, double outputVoltageMagnitude, double diodeVoltageDrop, double efficiency, double turnsRatio) {
+    double Cuk::calculate_duty_cycle(double inputVoltage, double outputVoltageMagnitude, double diodeVoltageDrop, double efficiency, double turnsRatio, double maximumDutyCycle) {
         // Cuk CCM, ideal:
         //   V1/V2 (n=1):     M(D) = -D/(1-D)
         //   V3 isolated:     M(D) = -D / ((1-D) · n),  n = Np/Ns (Flyback convention)
@@ -37,10 +37,14 @@ namespace OpenMagnetics {
         double effVin = inputVoltage * efficiency;
         double reflectedVo = (outputVoltageMagnitude + diodeVoltageDrop) * turnsRatio;
         double dutyCycle = reflectedVo / (reflectedVo + effVin);
-        if (dutyCycle >= 0.95) {
+        // Explicit configurable maxD gate (no silent clamp). 1 % tolerance
+        // absorbs design-requirements rounding.
+        constexpr double dutyTolerance = 0.01;
+        if (dutyCycle > maximumDutyCycle * (1.0 + dutyTolerance)) {
             throw InvalidInputException(ErrorCode::INVALID_INPUT,
                 "Cuk::calculate_duty_cycle: duty cycle " + std::to_string(dutyCycle) +
-                " >= 0.95 — converter would lose regulation; reduce |Vo| or raise Vin");
+                " exceeds maximumDutyCycle " + std::to_string(maximumDutyCycle) +
+                " — converter would lose regulation; reduce |Vo|, raise Vin, or relax maximumDutyCycle.");
         }
         return dutyCycle;
     }
@@ -114,7 +118,7 @@ namespace OpenMagnetics {
             for (const auto& op : get_operating_points()) {
                 double Iout = op.get_output_currents()[0];
                 double Vo   = std::abs(op.get_output_voltages()[0]);
-                double D    = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio);
+                double D    = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio, maximumDutyCycle.value_or(0.95));
                 double IL1avg = Iout * D / ((1.0 - D) * turnsRatio * efficiency);
                 maximumDeltaIL1 = std::max(maximumDeltaIL1, rippleRatio * IL1avg);
             }
@@ -124,7 +128,7 @@ namespace OpenMagnetics {
             for (const auto& op : get_operating_points()) {
                 double Iout = op.get_output_currents()[0];
                 double Vo   = std::abs(op.get_output_voltages()[0]);
-                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio);
+                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio, maximumDutyCycle.value_or(0.95));
                 double IL1avg = Iout * D / ((1.0 - D) * turnsRatio * efficiency);
                 double IL2avg = Iout;
                 // Switch-current peak: IS_pk = IL1avg + IL2avg/n + ripples/2.
@@ -141,7 +145,7 @@ namespace OpenMagnetics {
         for (const auto& op : get_operating_points()) {
             double switchingFrequency = op.get_switching_frequency();
             double Vo = std::abs(op.get_output_voltages()[0]);
-            double D  = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio);
+            double D  = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio, maximumDutyCycle.value_or(0.95));
             double L1 = calculate_l1_min(maximumInputVoltage, D, maximumDeltaIL1, switchingFrequency);
             maximumNeededInductance = std::max(maximumNeededInductance, L1);
         }
@@ -186,7 +190,7 @@ namespace OpenMagnetics {
         }
         lastTurnsRatio = turnsRatio;
 
-        double dutyCycle = calculate_duty_cycle(inputVoltage, outputVoltageMag, diodeVoltageDrop, efficiency, turnsRatio);
+        double dutyCycle = calculate_duty_cycle(inputVoltage, outputVoltageMag, diodeVoltageDrop, efficiency, turnsRatio, maximumDutyCycle.value_or(0.95));
 
         // Common derived quantities (V1/V2/V3).
         // Power balance: Vin·η·IL1avg = |Vo|·Iout  ⇒  IL1avg = |Vo|·Iout/(Vin·η).
@@ -467,7 +471,7 @@ namespace OpenMagnetics {
             for (const auto& op : get_operating_points()) {
                 double Iout = op.get_output_currents()[0];
                 double Vo   = std::abs(op.get_output_voltages()[0]);
-                double D    = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency);
+                double D    = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, 1.0, maximumDutyCycle.value_or(0.95));
                 double IL1avg = Iout * D / (1.0 - D);
                 maximumDeltaIL1 = std::max(maximumDeltaIL1, rippleRatio * IL1avg);
             }
@@ -477,7 +481,7 @@ namespace OpenMagnetics {
             for (const auto& op : get_operating_points()) {
                 double Iout = op.get_output_currents()[0];
                 double Vo   = std::abs(op.get_output_voltages()[0]);
-                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency);
+                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, 1.0, maximumDutyCycle.value_or(0.95));
                 double IL1avg = Iout * D / (1.0 - D);
                 // IS_peak ≈ IL1avg + IL2avg + (ΔIL1 + ΔIL2)/2; treat ΔIL2 as
                 // ~30 % of IL2avg (the L2 sizing default), so the L1 ripple
@@ -496,7 +500,7 @@ namespace OpenMagnetics {
         for (const auto& op : get_operating_points()) {
             double switchingFrequency = op.get_switching_frequency();
             double Vo = std::abs(op.get_output_voltages()[0]);
-            double D  = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency);
+            double D  = calculate_duty_cycle(maximumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, 1.0, maximumDutyCycle.value_or(0.95));
             double L1 = calculate_l1_min(maximumInputVoltage, D, maximumDeltaIL1, switchingFrequency);
             maximumNeededInductance = std::max(maximumNeededInductance, L1);
         }
@@ -538,7 +542,7 @@ namespace OpenMagnetics {
                 double switchingFrequency = op.get_switching_frequency();
                 double Vo   = std::abs(op.get_output_voltages()[0]);
                 double Iout = op.get_output_currents()[0];
-                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio);
+                double D    = calculate_duty_cycle(minimumInputVoltage, Vo, get_diode_voltage_drop(), efficiency, turnsRatio, maximumDutyCycle.value_or(0.95));
                 double IL1avg = Iout * D / ((1.0 - D) * turnsRatio * efficiency);
                 double VCa    = minimumInputVoltage / (1.0 - D);
                 double deltaIm_target = std::max(0.20 * IL1avg, 1e-6);
@@ -719,7 +723,7 @@ namespace OpenMagnetics {
         double efficiency = 1.0;
         if (get_efficiency()) efficiency = get_efficiency().value();
 
-        double dutyCycle = calculate_duty_cycle(inputVoltage, outputVoltageMag, diodeVoltageDrop, efficiency);
+        double dutyCycle = calculate_duty_cycle(inputVoltage, outputVoltageMag, diodeVoltageDrop, efficiency, 1.0, maximumDutyCycle.value_or(0.95));
 
         // Internally-sized L2, C1, Co (mirror process_operating_points_for_input_voltage)
         double IL1avg = outputCurrent * dutyCycle / (1.0 - dutyCycle);

src/converter_models/Cuk.h

@@ -134,6 +134,19 @@ class Cuk : public MAS::Cuk, public Topology {
     int numPeriodsToExtract = 5;
     int numSteadyStatePeriods = 50;
 
+    // Maximum operating duty cycle. Cuk CCM duty is
+    // D = (|Vo|+Vd)*n / (Vin*eta + (|Vo|+Vd)*n); for low Vin or high
+    // step-down/up ratio this approaches 1, where the conversion ratio
+    // M = -D/(1-D) blows up and regulation collapses. Cuk previously
+    // hard-coded a 0.95 ceiling in calculate_duty_cycle(); this field
+    // makes it configurable so callers with stricter controller
+    // constraints (e.g. maxD = 0.7) get a loud throw rather than
+    // running well past their controller's range. Mirrors Flyback
+    // (04272d7b), forward family (683e731c), Buck (2c9300c2), Boost
+    // (96fdb52a), IsolatedBuck (703bc80e), IsolatedBuckBoost
+    // (a158d548).
+    std::optional<double> maximumDutyCycle = 0.95;
+
     // Internal sizing rules-of-thumb for L2 / C1 / Co (V1 only).
     // These are fixed defaults; future revisions will expose them via
     // schema fields and/or DesignRequirements (see CUK_PLAN.md §13).
@@ -199,6 +212,9 @@ class Cuk : public MAS::Cuk, public Topology {
     int get_num_steady_state_periods() const { return numSteadyStatePeriods; }
     void set_num_steady_state_periods(int value) { this->numSteadyStatePeriods = value; }
 
+    std::optional<double> get_maximum_duty_cycle() const { return maximumDutyCycle; }
+    void set_maximum_duty_cycle(std::optional<double> value) { this->maximumDutyCycle = value; }
+
     // ---- Per-OP diagnostic accessors ----
     double get_last_duty_cycle()                  const { return lastDutyCycle; }
     double get_last_conversion_ratio()            const { return lastConversionRatio; }
@@ -251,7 +267,7 @@ class Cuk : public MAS::Cuk, public Topology {
     //
     // For V3 isolated, turnsRatio = Np/Ns (Flyback convention). Conversion gain
     // becomes |Vo|/Vin = D / ((1-D) · turnsRatio); reduces to V1 when turnsRatio=1.
-    static double calculate_duty_cycle(double inputVoltage, double outputVoltageMagnitude, double diodeVoltageDrop, double efficiency, double turnsRatio = 1.0);
+    static double calculate_duty_cycle(double inputVoltage, double outputVoltageMagnitude, double diodeVoltageDrop, double efficiency, double turnsRatio = 1.0, double maximumDutyCycle = 0.95);
     static double calculate_conversion_ratio(double dutyCycle);                 // -D/(1-D)
     static double calculate_coupling_cap_voltage(double inputVoltage, double dutyCycle);
     static double calculate_l1_min(double inputVoltage, double dutyCycle, double deltaIL1, double switchingFrequency);