# blast_load.py -- Kinney-Graham serbest alan patlama parametreleri + normal yansima
# Kaynak: Kinney & Graham (1985), Explosive Shocks in Air; Rankine-Hugoniot (gamma=1.4)
import math

W  = 100.0      # kg TNT esdegeri
R  = 20.0       # m, mesafe
P0 = 101.325    # kPa, atmosfer basinci (deniz seviyesi)

Wc = W ** (1/3)
Z  = R / Wc     # m/kg^(1/3), olcekli mesafe

# Tepe serbest alan asinc basinci
num = 808 * (1 + (Z/4.5)**2)
den = math.sqrt((1+(Z/0.048)**2) * (1+(Z/0.32)**2) * (1+(Z/1.35)**2))
Pso = P0 * num / den                                   # kPa

# Pozitif faz suresi -- DOGRU bicim: YALNIZCA son terim karekok icinde
tn  = 980 * (1 + (Z/0.54)**10)
tdn = (1+(Z/0.02)**3) * (1+(Z/0.74)**6) * math.sqrt(1+(Z/6.9)**2)
td  = tn / tdn * Wc                                    # ms

# HATALI bicim (eski v2.0 kodu): (Z/0.74)^6 terimi de karekok icine alinmis
tdn_bad = (1+(Z/0.02)**3) * math.sqrt((1+(Z/0.74)**6) * (1+(Z/6.9)**2))
td_bad  = tn / tdn_bad * Wc

# Normal yansima (Rankine-Hugoniot, ideal gaz)
Cr = 2 * (7*P0 + 4*Pso) / (7*P0 + Pso)
Pr = Cr * Pso

print(f"Z        = {Z:.4f} m/kg^(1/3)")
print(f"Pso      = {Pso:.2f} kPa")
print(f"td       = {td:.3f} ms   (dogru bicim)")
print(f"td_hatali= {td_bad:.3f} ms   (eski kod) -> {td_bad/td:.1f} kat hata")
print(f"Cr       = {Cr:.4f}")
print(f"Pr       = {Pr:.2f} kPa")