//==================================================================================================
/*
  ROTGEN - Runtime Overlay for Eigen
  Copyright : CODE RECKONS
  SPDX-License-Identifier: BSL-1.0
*/
//==================================================================================================
#include <rotgen/rotgen.hpp>

#include "unit/tests.hpp"

TTS_CASE_TPL("SVD decomposition - Dynamic case", rotgen::tests::types)

<typename T, typename O>(tts::type<tts::types<T, O>>)
{
  int rank, i = 5;
  auto eps = std::numeric_limits<T>::epsilon();

  auto m = rotgen::setRandom<
    rotgen::matrix<T, rotgen::Dynamic, rotgen::Dynamic, O::value>>(5, 5);
  auto decomp = rotgen::svd(m);

  do
  {
    rank = decomp.rank();

    auto u = decomp.matrixU(rank);
    auto d = decomp.singularValues(rank);
    auto dd = decomp.matrixD(rank);
    auto v = decomp.matrixV(rank);

    TTS_EQUAL(rank, i);

    auto rec = (u * dd * rotgen::transpose(v));
    auto error = m - rec;

    TTS_LESS_EQUAL(rotgen::maxCoeff(rotgen::abs(error)) / eps, 16.)
      << "Result:\n"
      << rec << "\n"
      << "Residuals:\n"
      << error << "\n";

    // Reduce rank by duplicating one column
    i--;
    col(m, i) = col(m, 0);
    decomp = rotgen::svd(m);

  } while (rank != 1);
};

TTS_CASE_TPL("SVD decomposition - Static case", rotgen::tests::types)

<typename T, typename O>(tts::type<tts::types<T, O>>)
{
  int rank, i = 5;
  auto eps = std::numeric_limits<T>::epsilon();

  auto m = rotgen::setRandom<rotgen::matrix<T, 5, 5, O::value>>();
  auto decomp = rotgen::svd(m);

  do
  {
    rank = decomp.rank();

    auto u = decomp.matrixU(rank);
    auto d = decomp.singularValues(rank);
    auto dd = decomp.matrixD(rank);
    auto v = decomp.matrixV(rank);

    TTS_EQUAL(rank, i);

    auto rec = (u * dd * rotgen::transpose(v));
    auto error = m - rec;

    TTS_LESS_EQUAL(rotgen::maxCoeff(rotgen::abs(error)) / eps, 16.)
      << "Result:\n"
      << rec << "\n"
      << "Residuals:\n"
      << error << "\n";

    // Reduce rank by duplicating one column
    i--;
    col(m, i) = col(m, 0);
    decomp = rotgen::svd(m);

  } while (rank != 1);
};