#include "fiko/dokum/gerilme.hpp"

#include <algorithm>
#include <cmath>

namespace fiko::dokum {

namespace {

bool metalMi(const VoxelKafes& k, u32 idx) {
  return k.tip[idx] == HucreTip::Ic || k.tip[idx] == HucreTip::Sinir;
}

// G(fs) = 1/(1−fs) + 2·ln(1−fs) − (1−fs) : ∫ fs²/(1−fs)² dfs ilkeli.
f64 rdgIlkel(f64 fs) {
  const f64 u = 1.0 - fs;
  if (u <= 1e-9) return 0.0;  // fs→1 tekilliği — pencere üstü 0.99'la sınırlı
  return 1.0 / u + 2.0 * std::log(u) - u;
}

}  // namespace

f64 termalGenlesme(const Alasim& alasim) {
  if (alasim.tLikidus <= 800.0) return 23e-6;   // Al alaşımları
  if (alasim.tLikidus <= 1350.0) return 10e-6;  // gri dökme demir (GG)
  return 12e-6;                                 // döküm çeliği
}

f64 sicaklikFs(const Alasim& alasim, f64 fsHedef) {
  if (fsHedef <= 0.0) return alasim.tLikidus;
  if (fsHedef >= 1.0) return alasim.tSolidus;
  f64 alt = alasim.tSolidus, ust = alasim.tLikidus;  // fs(alt)=1, fs(ust)=0
  for (int i = 0; i < 64; ++i) {
    const f64 orta = 0.5 * (alt + ust);
    if (alasim.katiOran(orta) > fsHedef) alt = orta;  // fs T'de azalan
    else ust = orta;
  }
  return 0.5 * (alt + ust);
}

f64 cscAlasimOrani(const Alasim& alasim, const GerilmeAyar& ayar) {
  const f64 dTv = sicaklikFs(alasim, ayar.fsHassasAlt) -
                  sicaklikFs(alasim, ayar.fsHassasUst);
  const f64 dTr = sicaklikFs(alasim, ayar.fsRahatAlt) -
                  sicaklikFs(alasim, ayar.fsRahatUst);
  if (!(dTv > 0.0) || !(dTr > 1e-12)) return 0.0;
  return dTv / dTr;
}

f64 rdgDrenajFaktoru(const GerilmeAyar& ayar) {
  if (!(ayar.fsHassasUst > ayar.fsHassasAlt) || ayar.fsHassasUst >= 1.0)
    return 0.0;
  return rdgIlkel(ayar.fsHassasUst) - rdgIlkel(ayar.fsHassasAlt);
}

GerilmeSonuc gerilmeHesapla(const VoxelKafes& kafes, const KatilasmaSonuc& kat,
                            const Alasim& alasim, const GerilmeAyar& ayar) {
  GerilmeSonuc g;
  if (!kafes.gecerliMi()) return g;
  const u32 N = kafes.hucreSayi();
  if (!kat.hazir || kat.katilasmaSuresi.size() != N || kat.sdas.size() != N ||
      kat.niyama.size() != N || kafes.tip.size() != N)
    return g;
  const f64 dTf = alasim.tLikidus - alasim.tSolidus;
  if (!(dTf > 1e-9)) return g;  // ötektik (donma aralığı yok) — pencere kurulamaz

  const f64 dTv = sicaklikFs(alasim, ayar.fsHassasAlt) -
                  sicaklikFs(alasim, ayar.fsHassasUst);
  const f64 kappa = cscAlasimOrani(alasim, ayar);
  const f64 F = rdgDrenajFaktoru(ayar);
  const f64 alfa = termalGenlesme(alasim);
  const f64 dTcek = std::max(alasim.tSolidus - ayar.ortamSicaklik, 0.0);

  g.csc.assign(N, kKatmanGecersiz);
  g.rdg.assign(N, kKatmanGecersiz);
  g.carpilma.assign(N, kKatmanGecersiz);

  const u32 nx = kafes.nx(), ny = kafes.ny(), nz = kafes.nz();
  const auto katVakti = [&](u32 idx) -> f64 {
    return (metalMi(kafes, idx)) ? kat.katilasmaSuresi[idx] : -1.0;
  };

  // 1. geçiş: t_maks + hassas-pencere süresi t_v + komşu gecikme farkı Δt.
  f64 tsMaks = 0.0;
  for (u32 i = 0; i < N; ++i) tsMaks = std::max(tsMaks, katVakti(i));

  std::vector<f64> tv(N, -1.0), gecikme(N, -1.0);
  std::vector<f64> tvGecerli;
  u32 iMin = nx, iMaks = 0, jMin = ny, jMaks = 0, kMin = nz, kMaks = 0;
  const int di[6] = {-1, 1, 0, 0, 0, 0};
  const int dj[6] = {0, 0, -1, 1, 0, 0};
  const int dk[6] = {0, 0, 0, 0, -1, 1};
  for (u32 kk = 0; kk < nz; ++kk)
    for (u32 j = 0; j < ny; ++j)
      for (u32 i = 0; i < nx; ++i) {
        const u32 idx = kafes.indeks(i, j, kk);
        const f64 ts = katVakti(idx);
        if (ts < 0.0) continue;
        iMin = std::min(iMin, i); iMaks = std::max(iMaks, i);
        jMin = std::min(jMin, j); jMaks = std::max(jMaks, j);
        kMin = std::min(kMin, kk); kMaks = std::max(kMaks, kk);
        // t_v = ΔT_hassas/Ṙ (sabit-Ṙ yeniden-kurulum; Ṙ SDAS'tan geri kazanım).
        const f64 R = sogumaHiziSdas(kat.sdas[idx], alasim);
        if (R > 1e-12 && dTv > 0.0) {
          tv[idx] = dTv / R;
          tvGecerli.push_back(tv[idx]);
        }
        // Komşu katılaşma-zamanı uyumsuzluğu (çarpılma/RDG gerinim vekili).
        f64 dt = 0.0;
        for (int f = 0; f < 6; ++f) {
          const int ni = (int)i + di[f], nj = (int)j + dj[f], nk = (int)kk + dk[f];
          if (ni < 0 || nj < 0 || nk < 0 || ni >= (int)nx || nj >= (int)ny ||
              nk >= (int)nz)
            continue;
          const f64 tsK = katVakti(kafes.indeks((u32)ni, (u32)nj, (u32)nk));
          if (tsK >= 0.0) dt = std::max(dt, std::abs(ts - tsK));
        }
        gecikme[idx] = (tsMaks > 1e-12) ? dt / tsMaks : 0.0;
      }

  // Parça medyanı t_v (CSC uzamsallaştırma tabanı — banner'da belgeli).
  f64 tvMed = 0.0;
  if (!tvGecerli.empty()) {
    const size_t p = (tvGecerli.size() - 1) / 2;
    std::nth_element(tvGecerli.begin(), tvGecerli.begin() + p, tvGecerli.end());
    tvMed = tvGecerli[p];
  }

  // 2. geçiş: alanlar.
  f64 maksGecikme = 0.0;
  for (u32 idx = 0; idx < N; ++idx) {
    if (gecikme[idx] < 0.0) continue;  // metal değil / katılaşmamış
    maksGecikme = std::max(maksGecikme, gecikme[idx]);
    if (tv[idx] > 0.0 && tvMed > 1e-12)
      g.csc[idx] = (f32)(kappa * tv[idx] / tvMed);
    // RDG-tabanlı HCS: yalnız Niyama örneklenmiş hücre (sdas>0 ⇔ örneklendi).
    if (kat.sdas[idx] > 0.0) {
      const f64 epsMek = alfa * dTf * gecikme[idx];
      const f64 drenaj = std::max(kat.niyama[idx] / ayar.nyRef, ayar.nyMinOran);
      g.rdg[idx] = (f32)(epsMek * F / drenaj);
    }
    g.carpilma[idx] = (f32)(alfa * dTcek * gecikme[idx]);
  }

  // Parça geneli maks sapma (ilk-mertebe, kaba): κ_eğrilik = Δε/t → κ·L²/8.
  if (iMaks >= iMin) {
    const f64 hMm = kafes.h() * 1000.0;
    const f64 ex = (iMaks - iMin + 1) * hMm, ey = (jMaks - jMin + 1) * hMm,
              ez = (kMaks - kMin + 1) * hMm;
    const f64 L = std::max({ex, ey, ez}), t = std::min({ex, ey, ez});
    const f64 epsFark = alfa * dTcek * maksGecikme;
    if (t > 1e-9) g.maksSapmaMm = epsFark * L * L / (8.0 * t);
  }

  g.hazir = true;
  return g;
}

}  // namespace fiko::dokum