#include <gtest/gtest.h>

#include <array>
#include <cmath>
#include <limits>
#include <string>

#include "fesa/elements/thermal_shell.h"
#include "fesa/model/domain.h"

namespace {
using fesa::Vector3;
using Patch = std::array<Vector3, 4>;
const Patch kSquare{{{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}}};

fesa::Result<fesa::ThermalShell> Kernel(const Patch& patch) {
  return fesa::ThermalShell::Create({"part", 27, "0027"}, {{0, 1, 2, 3}}, patch,
                                    1., 1., 1., 1., {"mesh.inp", 17});
}

fesa::ModelDefinition Mesh(
    const std::vector<Vector3>& nodes,
    const std::vector<std::array<fesa::EntityIndex, 4>>& elements) {
  fesa::ModelDefinition model{};
  for (std::size_t i = 0; i < nodes.size(); ++i) {
    model.nodes.push_back(
        {{"part", static_cast<std::int64_t>(i + 1), std::to_string(i + 1)},
         nodes[i].Components(),
         {"mesh.inp", i + 1}});
  }
  model.thermal_materials.push_back({{"part", 1, "1"}, "heat", 1, 1, 1, {}});
  model.thermal_sections.push_back({{"part", 1, "1"}, "section", 1, 0, {}});
  for (std::size_t i = 0; i < elements.size(); ++i) {
    model.thermal_elements.push_back(
        {{"part", static_cast<std::int64_t>(i + 1), std::to_string(i + 1)},
         elements[i],
         0,
         {"mesh.inp", i + 20}});
  }
  model.thermal_step = fesa::ThermalStepDefinition{"heat", .1, 1., 10};
  return model;
}

bool Accepted(fesa::ModelDefinition model) {
  return fesa::Domain::Create(std::move(model)).HasValue();
}
}  // namespace

TEST(ThermalShellGeometry, CertifiesAffineDistortedSaddleAndSlenderPatches) {
  EXPECT_TRUE(Kernel(kSquare).HasValue());
  EXPECT_TRUE(
      Kernel({{{0, 0, 0}, {2, 0, 0}, {3, 1, 0}, {1, 1, 0}}}).HasValue());
  EXPECT_TRUE(
      Kernel({{{0, 0, 0}, {2, 0, 0}, {1.5, 1, 0}, {0, 1, 0}}}).HasValue());
  EXPECT_TRUE(Kernel({{{0, 0, 0},
                       {1, 0, 0},
                       {1, std::ldexp(1., -40), 0},
                       {0, std::ldexp(1., -40), 0}}})
                  .HasValue());
  for (double a : {.1, .5, 2.}) {
    EXPECT_TRUE(
        Kernel({{{-1, -1, a}, {1, -1, -a}, {1, 1, a}, {-1, 1, -a}}}).HasValue())
        << a;
  }
}

TEST(ThermalShellGeometry, PowerOfTwoScalingTranslationAndRotation) {
  for (int exponent : {-200, 0, 200}) {
    Patch patch = kSquare;
    for (auto& point : patch) point = std::ldexp(1., exponent) * point;
    EXPECT_TRUE(Kernel(patch).HasValue()) << exponent;
  }
  Patch translated = kSquare;
  Patch rotated = kSquare;
  for (std::size_t i = 0; i < 4; ++i) {
    translated[i] = translated[i] + Vector3{1048576., -1048576., 1048576.};
    rotated[i] = Vector3{.6 * kSquare[i].X() - .8 * kSquare[i].Y(),
                         .8 * kSquare[i].X() + .6 * kSquare[i].Y(), 0};
  }
  EXPECT_TRUE(Kernel(translated).HasValue());
  EXPECT_TRUE(Kernel(rotated).HasValue());
}

TEST(ThermalShellGeometry,
     RejectsWholePatchSingularityDespiteNonzeroGaussPoints) {
  Patch patch{};
  const std::array<double, 4> xi{{-1, 1, 1, -1}};
  const std::array<double, 4> eta{{-1, -1, 1, 1}};
  for (std::size_t i = 0; i < 4; ++i)
    patch[i] = {xi[i], (xi[i] - .25) * eta[i], 0};
  const auto result = Kernel(patch);
  ASSERT_FALSE(result.HasValue());
  EXPECT_EQ(result.GetStatus().Diagnostics()[0].entity_identity, "part:0027");
  EXPECT_EQ(result.GetStatus().Diagnostics()[0].location.line, 17U);
  for (double q : {-std::sqrt(.6), 0., std::sqrt(.6)})
    EXPECT_GT(std::abs(q - .25), 0.);
}

TEST(ThermalShellGeometry, RejectsInvalidAndDoesNotRepairUncertainGeometry) {
  EXPECT_FALSE(
      Kernel({{{0, 0, 0}, {1, 1, 0}, {1, 0, 0}, {0, 1, 0}}}).HasValue());
  EXPECT_FALSE(
      Kernel({{{0, 0, 0}, {1, 0, 0}, {2, 0, 0}, {3, 0, 0}}}).HasValue());
  Patch repeated = kSquare;
  repeated[2] = repeated[1];
  EXPECT_FALSE(Kernel(repeated).HasValue());
  for (double bad : {std::numeric_limits<double>::infinity(),
                     std::numeric_limits<double>::quiet_NaN()}) {
    Patch patch = kSquare;
    patch[0] = {bad, 0, 0};
    EXPECT_FALSE(Kernel(patch).HasValue());
  }
  const auto duplicate = fesa::ThermalShell::Create(
      {"part", 1, "1"}, {{0, 1, 1, 3}}, kSquare, 1, 1, 1, 1);
  EXPECT_FALSE(duplicate.HasValue());
  const double e = std::ldexp(1., -52);
  const auto uncertain =
      Kernel({{{-2, -2 - e, 0}, {0, -e, 0}, {2, 2 + e, 0}, {0, e, 0}}});
  if (!uncertain.HasValue())
    EXPECT_EQ(uncertain.GetStatus().Diagnostics()[0].code,
              "uncertain-thermal-shell-geometry");
}

TEST(ThermalShellGeometry, DomainEnforcesOppositeEdgesAndOneVertexLink) {
  const std::vector<Vector3> nodes{{0, 0, 0}, {1, 0, 0}, {1, 1, 0},
                                   {0, 1, 0}, {2, 0, 0}, {2, 1, 0}};
  EXPECT_TRUE(Accepted(Mesh(nodes, {{{0, 1, 2, 3}}, {{1, 4, 5, 2}}})));
  EXPECT_FALSE(Accepted(Mesh(nodes, {{{0, 1, 2, 3}}, {{2, 5, 4, 1}}})));
  EXPECT_FALSE(Accepted(Mesh(nodes, {{{0, 1, 2, 3}}, {{0, 1, 2, 3}}})));
  EXPECT_FALSE(Accepted(Mesh({{0, 0, 0},
                              {1, 0, 0},
                              {1, 1, 0},
                              {0, 1, 0},
                              {-1, -1, 0},
                              {0, -1, 0},
                              {-1, 0, 0}},
                             {{{0, 1, 2, 3}}, {{0, 6, 4, 5}}})));
  EXPECT_FALSE(
      Accepted(Mesh({{0, 0, 0},
                     {1, 0, 0},
                     {1, 1, 0},
                     {0, 1, 0},
                     {2, 0, 0},
                     {2, 1, 0},
                     {1, 0, 1},
                     {1, 1, 1}},
                    {{{0, 1, 2, 3}}, {{1, 4, 5, 2}}, {{1, 6, 7, 2}}})));
}

TEST(ThermalShellGeometry, DomainAllowsFoldAndDisconnectedCoordinateContact) {
  auto folded =
      Mesh({{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}, {1, 0, 1}, {1, 1, 1}},
           {{{0, 1, 2, 3}}, {{1, 4, 5, 2}}});
  folded.thermal_materials.push_back({{"part", 2, "2"}, "heat2", 10, 1, 1, {}});
  folded.thermal_sections.push_back({{"part", 2, "2"}, "section2", .5, 1, {}});
  folded.thermal_elements[1].section_index = 1;
  EXPECT_TRUE(Accepted(std::move(folded)));
  const std::vector<Vector3> nodes{{0, 0, 0},  {2, 0, 0}, {2, 1, 0},
                                   {0, 1, 0},  {0, 0, 0}, {1, 0, 0},
                                   {1, -1, 0}, {0, -1, 0}};
  auto disconnected = Mesh(nodes, {{{0, 1, 2, 3}}, {{4, 5, 6, 7}}});
  EXPECT_TRUE(Accepted(disconnected));
  for (std::size_t i = 4; i < 8; ++i)
    disconnected.nodes[i].source_id.instance_name = "other";
  disconnected.thermal_elements[1].source_id.instance_name = "other";
  EXPECT_TRUE(Accepted(disconnected));
}

TEST(ThermalShellGeometry,
     DomainRejectsSourceSharingPartialEdgeWithoutGlobalScan) {
  // The patches share the full 3->0 / 0->3 edge: connected vertex links remain
  // paths, while the additional 0->4 edge lies inside coarse 0->1.
  const std::vector<Vector3> nodes{{0, 0, 0}, {2, 0, 0}, {2, 1, 0},
                                   {0, 1, 0}, {1, 0, 0}, {1, 1, 1}};
  const auto model = Mesh(nodes, {{{0, 1, 2, 3}}, {{0, 3, 5, 4}}});
  const auto result = fesa::Domain::Create(model);
  ASSERT_FALSE(result.HasValue());
  EXPECT_EQ(result.GetStatus().Diagnostics()[0].code,
            "uncertain-thermal-shell-connection");
  auto invalid =
      Mesh({{0, 0, 0}, {1, 1, 0}, {1, 0, 0}, {0, 1, 0}}, {{{0, 1, 2, 3}}});
  EXPECT_FALSE(Accepted(invalid));
}

TEST(ThermalShellGeometry, SearchesSeparatingDirectionsBeyondCenterAndCorners) {
  // Exact raw certificate s=(0,-3,6) has min(s.w)=9; center has zero margin.
  EXPECT_TRUE(
      Kernel({{{3, -9, -3}, {1, 3, -3}, {-1, 3, 3}, {-3, 3, 3}}}).HasValue());
  // Exact s=(1.5,2.5,1) gives min(s.w)=3.5, whereas center/corner candidates
  // cannot separate the normal image. Edge/interior minimization is needed.
  EXPECT_TRUE(
      Kernel({{{-1, 4, -4}, {-1, -2, 2}, {3, 2, 4}, {-1, -4, -2}}}).HasValue());
}

TEST(ThermalShellGeometry,
     RejectsAliasedSourceIdentityAndAcceptsClosedVertexLink) {
  auto aliases =
      Mesh({{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}}, {{{0, 1, 2, 3}}});
  aliases.nodes[2].source_id = aliases.nodes[1].source_id;
  EXPECT_FALSE(Accepted(aliases));
  EXPECT_TRUE(Accepted(
      Mesh({{-1, -1, 0},
            {0, -1, 0},
            {1, -1, 0},
            {-1, 0, 0},
            {0, 0, 0},
            {1, 0, 0},
            {-1, 1, 0},
            {0, 1, 0},
            {1, 1, 0}},
           {{{0, 1, 4, 3}}, {{1, 2, 5, 4}}, {{3, 4, 7, 6}}, {{4, 5, 8, 7}}})));
}