refactor: 提取 CloudExpansionModel 独立模块

- RadiusAt(t)/DensityAt(t)/TimeToReach(r) 封装在 CloudExpansionModel
- ComputeEffectiveRadius 简化为委托调用
- GaussianPuffDispersion 重复逻辑待后续统一
This commit is contained in:
tian 2026-06-13 18:14:42 +08:00
parent 9c153a040b
commit 510403a66c
2 changed files with 68 additions and 43 deletions

View File

@ -0,0 +1,61 @@
using System;
using CounterDrone.Core.Models;
namespace CounterDrone.Core.Algorithms
{
/// <summary>云团膨胀模型——根据弹药和环境参数计算各阶段半径、密度、膨胀时间</summary>
public class CloudExpansionModel
{
private readonly AmmunitionSpec _ammo;
private readonly CombatScene _env;
public CloudExpansionModel(AmmunitionSpec ammo, CombatScene env)
{
_ammo = ammo;
_env = env;
}
/// <summary>Phase 1 爆轰初始半径 (m)</summary>
public float InitialRadius => 3.3f * (float)Math.Pow(Math.Max(0.01, (float)_ammo.BurstChargeKg), 0.32f);
public float TurbulentExpansionK => (float)_ammo.TurbulentExpansionK;
/// <summary>指定时刻湍流膨胀半径 R(t) = R₀ + k√t</summary>
public float RadiusAt(float elapsedSeconds)
{
if (elapsedSeconds <= 0) return InitialRadius;
bool inPhase3 = elapsedSeconds > 30f;
if (!inPhase3)
return InitialRadius + TurbulentExpansionK * (float)Math.Sqrt(elapsedSeconds);
// Phase 3: 高斯扩散
float windSpeed = (float)_env.WindSpeed;
float x = Math.Max(1f, windSpeed * (elapsedSeconds - 30f));
var cls = Kinematics.GetStabilityClass((WeatherType)_env.WeatherType, windSpeed);
float sY = Kinematics.SigmaY(cls, x);
float sZ = Kinematics.SigmaZ(cls, x);
float peakC = Kinematics.GaussianPeakConcentration((float)_ammo.SourceStrength, sY, sZ);
float effConc = (float)_ammo.EffectiveConcentration;
if (peakC <= effConc) return InitialRadius + TurbulentExpansionK * (float)Math.Sqrt(30f);
float sigma = (sY + sZ) / 2f;
float rPhase2 = InitialRadius + TurbulentExpansionK * (float)Math.Sqrt(30f);
return rPhase2 + sigma * (float)Math.Sqrt(2f * Math.Log(peakC / effConc));
}
/// <summary>指定时刻中心浓度</summary>
public float DensityAt(float elapsedSeconds)
{
float r = RadiusAt(elapsedSeconds);
float volume = (4f / 3f) * (float)Math.PI * r * r * r;
return volume > 0.001f ? (float)_ammo.SourceStrength / volume : (float)_ammo.CoreDensity;
}
/// <summary>达到指定半径所需时间 (s)</summary>
public float TimeToReach(float radius)
{
float r0 = InitialRadius;
if (radius <= r0) return 0;
float k = TurbulentExpansionK;
return (radius - r0) * (radius - r0) / (k * k);
}
}
}

View File

@ -349,41 +349,13 @@ namespace CounterDrone.Core.Algorithms
private (float effectiveRadius, float expansionTime, float rPhase2) ComputeEffectiveRadius(
DroneGroup threat, AmmunitionSpec ammo, CombatScene env)
{
float totalFlightTime = threat.ArrivalTime * 2f;
float halfTime = totalFlightTime / 2f;
float windSpeed = (float)env.WindSpeed;
var weather = (WeatherType)env.WeatherType;
float avgSpeed = (float)threat.Target.TypicalSpeed / 3.6f;
// Phase 1 初始半径(基于爆发药 TNT 当量)
float r0 = 3.3f * (float)Math.Pow(Math.Max(0.01, (float)ammo.BurstChargeKg), 0.32);
float k = (float)ammo.TurbulentExpansionK;
float maxDur = (float)ammo.MaxDuration;
// Phase 2 湍流膨胀后的半径
float phase2Time = Math.Min(halfTime, 30f);
float turbulentRadius = r0 + k * (float)Math.Sqrt(Math.Max(0, phase2Time));
float expansionTime = (float)Math.Pow((turbulentRadius - r0) / Math.Max(k, 0.01f), 2);
float effectiveR = turbulentRadius;
// Phase 3 高斯扩散
if (halfTime > 30f)
{
float x = Math.Max(1f, windSpeed * (halfTime - 30f));
var cls = Kinematics.GetStabilityClass(weather, windSpeed);
float sY = Kinematics.SigmaY(cls, x);
float sZ = Kinematics.SigmaZ(cls, x);
float peakC = Kinematics.GaussianPeakConcentration((float)ammo.SourceStrength, sY, sZ);
float threshold = (float)ammo.EffectiveConcentration;
if (peakC > threshold)
{
float sigma = (sY + sZ) / 2f;
effectiveR = turbulentRadius + sigma * (float)Math.Sqrt(2f * Math.Log(peakC / threshold));
}
}
return (effectiveR, expansionTime, turbulentRadius);
var model = new CloudExpansionModel(ammo, env);
float tPhase2 = 30f;
float rPhase2 = model.RadiusAt(tPhase2);
float expansionTime = model.TimeToReach(rPhase2);
float halfTime = threat.ArrivalTime * 2f; // 总飞行时间的一半
float effectiveR = model.RadiusAt(halfTime);
return (effectiveR, expansionTime, rPhase2);
}
private int CalcRoundsNeeded(DroneGroup threat, AmmunitionSpec ammo,
@ -399,14 +371,6 @@ namespace CounterDrone.Core.Algorithms
return Math.Max(1,
(int)Math.Ceiling(requiredCoverage / spacing));
var aerosolType = (AerosolType)ammo.AerosolType;
float neededExposure = aerosolType == AerosolType.ActiveMaterial ? 2f : 6f;
float singleCoverage = neededExposure * avgSpeed;
float requiredCoverage = singleCoverage;
return Math.Max(1,
(int)Math.Ceiling((requiredCoverage - 2f * turbulentRadius) / spacing) + 1);
}
private float ComputeInterceptProbability(DroneGroup threat,