diff --git a/src/gausse/gpu/batch_integrator.py b/src/gausse/gpu/batch_integrator.py index 396a325..3c8baac 100644 --- a/src/gausse/gpu/batch_integrator.py +++ b/src/gausse/gpu/batch_integrator.py @@ -108,7 +108,8 @@ def _derivatives(xp, q, i, x, v, p: BatchDischargeParams): v_c = q / p.capacitance_f d_q = -i d_i = (v_c - i * r_eff - dlambda_dx * v) / dlambda_di - retard = xp.sign(v) * (p.retard_const_n + p.drag_coeff_n * v * v) # трение + воздух + # трение + воздух; tanh-сглаживание у v=0 — как в physics/circuit.py + retard = xp.tanh(v / 0.01) * (p.retard_const_n + p.drag_coeff_n * v * v) d_v = (force - retard) / p.mass_kg d_x = v return d_q, d_i, d_x, d_v @@ -145,7 +146,7 @@ def _get_fused_step(cp): reff = rt + red * ov dq = -i di = (q / C - i * reff - liron * dov * g * v) / dl_di - dv = (liron * dov * G - cp.sign(v) * (fc + dc * v * v)) / m + dv = (liron * dov * G - cp.tanh(v * 100.0) * (fc + dc * v * v)) / m return dq, di, v, dv @cp.fuse() diff --git a/src/gausse/physics/circuit.py b/src/gausse/physics/circuit.py index 2d96321..f7f6e93 100644 --- a/src/gausse/physics/circuit.py +++ b/src/gausse/physics/circuit.py @@ -35,9 +35,18 @@ class StageCircuitParams: drag_coeff_n_per_mps2: float = 0.0 +# Сглаживание сухого трения около v=0: sign(v) -> tanh(v/ε). РАЗРЫВНАЯ сила +# (чистый sign) ломает адаптивный решатель: у v≈0 шаг дробится бесконечно +# (chattering), одна конфигурация считается десятки минут. tanh — стандартная +# регуляризация (модель залипания); при |v| >> ε неотличима от sign. +FRICTION_SMOOTHING_V_MPS = 0.01 + + def retarding_force_n(params: StageCircuitParams, v: float) -> float: - """Тормозящая сила (трение + воздух) со знаком ПРОТИВ скорости.""" - return float(np.sign(v)) * (params.retard_const_n + params.drag_coeff_n_per_mps2 * v * v) + """Тормозящая сила (трение + воздух), гладкая по v, против скорости.""" + return float(np.tanh(v / FRICTION_SMOOTHING_V_MPS)) * ( + params.retard_const_n + params.drag_coeff_n_per_mps2 * v * v + ) def effective_resistance(inductance_model: CoilInductanceModel, params: StageCircuitParams, x) -> float: diff --git a/src/gausse/sim/stage.py b/src/gausse/sim/stage.py index 04d927e..e9dce72 100644 --- a/src/gausse/sim/stage.py +++ b/src/gausse/sim/stage.py @@ -122,10 +122,14 @@ def _switch_equivalent_resistance_ohm( def _coast_derivatives(t: float, state: np.ndarray, retard_const_n: float, drag_coeff: float, mass_kg: float) -> list[float]: - """Подлёт к датчику: не баллистика в вакууме, а с трением о трубку и воздухом.""" + """Подлёт к датчику: не баллистика в вакууме, а с трением о трубку и воздухом. + + Трение сглажено tanh'ом (см. circuit.FRICTION_SMOOTHING_V_MPS) — разрывный + sign(v) заставляет адаптивный решатель бесконечно дробить шаг у v≈0. + """ v = state[1] - retard = retard_const_n + drag_coeff * v * v - return [v, -math.copysign(retard, v) / mass_kg if v != 0.0 else 0.0] + retard = (retard_const_n + drag_coeff * v * v) * math.tanh(v / circuit.FRICTION_SMOOTHING_V_MPS) + return [v, -retard / mass_kg] def _stall_event(t: float, state: np.ndarray, *_args) -> float: diff --git a/tests/test_realism_v2.py b/tests/test_realism_v2.py index 153d473..17b7521 100644 --- a/tests/test_realism_v2.py +++ b/tests/test_realism_v2.py @@ -145,3 +145,22 @@ def test_pulse_limit_fallback_without_datasheet(): sw = replace(DB.switches[0], pulse_current_a=None) # без даташита — консервативный множитель по типу ключа assert switch_pulse_limit_a(sw) == sw.max_current_a * {"SCR": 10.0, "MOSFET": 4.0, "IGBT": 3.0}[sw.kind] + + +def test_friction_force_smooth_at_zero_velocity(): + """Регрессия на завис решателя: сила трения обязана быть ГЛАДКОЙ в v=0. + + Разрывный sign(v) заставлял solve_ivp бесконечно дробить шаг у v≈0 + (одна конфигурация считалась 35+ минут на сервере). + """ + from gausse.physics.circuit import StageCircuitParams, retarding_force_n + + p = StageCircuitParams( + capacitance_f=1e-4, r_total_ohm=0.1, mass_kg=0.01, + retard_const_n=0.034, drag_coeff_n_per_mps2=2e-5, + ) + # непрерывность: около нуля сила ~ 0, нечётная, без скачка + assert abs(retarding_force_n(p, 1e-6)) < 1e-5 + assert abs(retarding_force_n(p, 1e-6) + retarding_force_n(p, -1e-6)) < 1e-12 + # вдали от нуля выходит на полное трение + assert retarding_force_n(p, 1.0) > 0.9 * p.retard_const_n