// GenFem - A high-level finite element library
// Copyright (C) 2010-2026 Eric Bechet
//
// See the LICENSE file for license information and contributions.
// Please report all bugs and problems to <bechet@cadxfem.org>.
//
// Initial design: Frederic Duboeuf and Eric Bechet (rev.576)

#ifndef _GEN_MATRICIAL_TRAITS__H_
#define _GEN_MATRICIAL_TRAITS__H_

#include <ostream>
#include <iomanip>

#include "genTensorialTraits.h"
#include "genTensorAdd.h"
#include "genMatrix.h"

//--------------------------------------------------------------------------
// Matricial Traits
//--------------------------------------------------------------------------

template<class T, int n> struct ContainerTraits
{
};


template<class T> struct ContainerTraits<T,0>
{
  typedef genScalar<T> Container;
};

template<class T> struct ContainerTraits<T,1>
{
  typedef genVector<T> Container;
};

template<class T> struct ContainerTraits<T,2>
{
  typedef genMatrix<T> Container;
};


//--------------------------------------------------------------------------
// Container Operators
//--------------------------------------------------------------------------

template<class T> void InitContainer(genScalar<T> &cont, int n1, int n2)
{
  cont=T();
}

template<class T> void InitContainer(genVector<T> &cont, int n1, int n2)
{
  cont.resize(n1);
}

template<class T> void InitContainer(genMatrix<T> &cont, int n1, int n2)
{
  cont.resize(n1,n2);
}


template<class T> void Print(const genScalar<T> &cont, const char* name="")
{
  printf("Container 0 : %s\n", name);
  PrintOp(cont());
}

template<class T> void Print(const genVector<T> &cont, const char* name="")
{
  printf("Container 1 : %s\n", name);
  int nbFF = cont.size();
  printf(" dim : %d\n", nbFF);
  for(int i = 0; i < nbFF; ++i)
  {
    printf(" FF(%d)\n", i);
    PrintOp(cont(i));
  }
}

template<class T> void Print(const genMatrix<T> &cont, const char* name="")
{
  printf("Container 2 : %s\n", name);
  int nbFF0 = cont.size1();
  int nbFF1 = cont.size2();
  printf(" dim : (%d,%d)\n", nbFF0, nbFF1);
  for(int i = 0; i < nbFF0; ++i)
    for(int j = 0; j < nbFF1; ++j)
    {
      printf(" FF(%d,%d)\n", i,j);
      PrintOp(cont(i,j));
    }
}


template<class T> void writeContainer(const genScalar<T> &cont, std::ostream &out)
{
  out << "ContainerDim 0 : \n";
  out << "\n";
  out.precision(5);
  out << std::setiosflags( std::ios::showpos );
  out << std::setiosflags( std::ios::scientific );
  int tensLength = TensorialProperties<T>::Length();
  for(int l = 0; l < tensLength; ++l)
  {
    out << "   ";
    writeTensor(cont(),out,l);
    out << "\n";
  }
  out << "\n";
  out << std::resetiosflags( std::ios::showpos );
  out << std::resetiosflags( std::ios::scientific );
}

template<class T> void writeContainer(const genVector<T> &cont, std::ostream &out)
{
  int nbFF = cont.size();
  out << "ContainerDim 1 : Size " << nbFF << "\n";
  out << "\n";
  out.precision(5);
  out << std::setiosflags( std::ios::showpos );
  out << std::setiosflags( std::ios::scientific );
  int tensLength = TensorialProperties<T>::Length();
  for(int i = 0; i < nbFF; ++i)
  {
    for(int l = 0; l < tensLength; ++l)
    {
      out << "   ";
      writeTensor(cont(i),out,l);
      out << "\n";
    }
    out << "\n";
  }
  out << "\n";
  out << std::resetiosflags( std::ios::showpos );
  out << std::resetiosflags( std::ios::scientific );
}

template<class T> void writeContainer(const genMatrix<T> &cont, std::ostream &out)
{
  int nbFF1 = cont.size1();
  int nbFF2 = cont.size2();
  out << "ContainerDim 2 : Size (" << nbFF1 << "," << nbFF2 << ") \n";
  out << "\n";
  out.precision(5);
  out << std::setiosflags( std::ios::showpos );
  out << std::setiosflags( std::ios::scientific );
  int tensLength = TensorialProperties<T>::Length();
  for(int i = 0; i < nbFF1; ++i)
  {
    for(int l = 0; l < tensLength; ++l)
    {
      out << "   ";
      for(int j = 0; j < nbFF2; ++j)
      {
        writeTensor(cont(i,j),out,l);
      }
      out << "\n";
    }
    out << "\n";
  }
  out << "\n";
  out << std::resetiosflags( std::ios::showpos );
  out << std::resetiosflags( std::ios::scientific );
}

template<class T> void SpaceUnion(const bool same[2], const genVector<T> &cont, const genVector<T> &val, genVector<T> &res)
{
  res.resize(cont.size()+val.size(),false);
  int i,s=cont.size();
  for (i=0; i<cont.size(); ++i) res(i)=cont(i);
  for (i=0; i<val.size(); ++i) res(i+s)=val(i);
//  std::cout << res.size() << "v" << std::endl;
}

template<class T> void SpaceUnion(const bool same[2], const genMatrix<T> &cont, const genMatrix<T> &val, genMatrix<T> &res)
{
  if ((same[1])&&(!same[0]))
  {
    res.resize(cont.size1()+val.size1(),cont.size2(),false);// no need to init the values to 0 because cont(i,j) and val(i,j) cover exactly res(i,j)
//    std::cout << res.size1() << " " << res.size2() << " m1" << std::endl;
    int i,j,s1=cont.size1();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i+s1,j)=val(i,j);
    return;
  }
  if ((!same[1])&&(same[0]))
  {
    res.resize(cont.size1(),cont.size2()+val.size2(),false);// no need to init the values to 0 because cont(i,j) and val(i,j) cover exactly res(i,j)
//        std::cout << res.size1() << " " << res.size2() << " m2" << std::endl;
    int i,j,s2=cont.size2();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i,j+s2)=val(i,j);
    return;
  }
  if ((!same[1])&&(!same[0]))
  {
    res.resize(cont.size1()+val.size1(),cont.size2()+val.size2()); // should init the values to 0 because cont(i,j) and val(i,j) don't cover res(i,j) completely
    // this is clearly suboptimal and here exclusively to let the program function even in that case. Maybe one should throw() ?
//        std::cout << res.size1() << " " << res.size2() << " m3" << std::endl;
    int i,j,s1=cont.size1(),s2=cont.size2();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i+s1,j+s2)=val(i,j);
    return;
  }
}



template<class T> void SpaceUnionm(const bool same[2], const genVector<T> &cont, const genVector<T> &val, genVector<T> &res)
{
  res.resize(cont.size()+val.size(),false);
  int i,s=cont.size();
  for (i=0; i<cont.size(); ++i) res(i)=cont(i);
  for (i=0; i<val.size(); ++i) res(i+s)=-val(i);
//  std::cout << res.size() << "v" << std::endl;
}

template<class T> void SpaceUnionm(const bool same[2], const genMatrix<T> &cont, const genMatrix<T> &val, genMatrix<T> &res)
{
  if ((same[1])&&(!same[0]))
  {
    res.resize(cont.size1()+val.size1(),cont.size2(),false);// no need to init the values to 0 because cont(i,j) and val(i,j) cover exactly res(i,j)
//    std::cout << res.size1() << " " << res.size2() << " m1" << std::endl;
    int i,j,s1=cont.size1();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i+s1,j)=-val(i,j);
    return;
  }
  if ((!same[1])&&(same[0]))
  {
    res.resize(cont.size1(),cont.size2()+val.size2(),false);// no need to init the values to 0 because cont(i,j) and val(i,j) cover exactly res(i,j)
//        std::cout << res.size1() << " " << res.size2() << " m2" << std::endl;
    int i,j,s2=cont.size2();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i,j+s2)=-val(i,j);
    return;
  }
  if ((!same[1])&&(!same[0]))
  {
    res.resize(cont.size1()+val.size1(),cont.size2()+val.size2()); // should init the values to 0 because cont(i,j) and val(i,j) don't cover res(i,j) completely
    // this is clearly suboptimal and here exclusively to let the program function even in that case. Maybe one should throw() ?
//        std::cout << res.size1() << " " << res.size2() << " m3" << std::endl;
    int i,j,s1=cont.size1(),s2=cont.size2();
    for (i=0; i<cont.size1(); ++i)
      for (j=0; j<cont.size2(); ++j)
        res(i,j)=cont(i,j);
    for (i=0; i<val.size1(); ++i)
      for (j=0; j<val.size2(); ++j)
        res(i+s1,j+s2)=-val(i,j);
    return;
  }
}





//--------------------------------------------------------------------------
// Shortcut Operators
//--------------------------------------------------------------------------

template<class T> void Cat(const genScalar<T> &cont, const genScalar<T> &val, genScalar<T> &res,const bool same[2])
{
   PlusOp<T>()(cont(),val(),res());
}

template<class T> void Cat(const genVector<T> &cont, const genVector<T> &val, genVector<T> &res,const bool same[2])
{
  if (same[0])
  {
    res.resize(cont.size(),false);
    for (int i=0;i<cont.size();++i)
      PlusOp<T>()(cont(i),val(i),res(i));
  }
  else SpaceUnion(same,cont,val,res);
}

template<class T> void Cat(const genMatrix<T> &cont, const genMatrix<T> &val, genMatrix<T> &res,const bool same[2])
{
  if ((same[0])&&(same[1]))
  {
    res.resize(cont.size1(), cont.size2(),false);
    for (int i=0;i<cont.size1();++i)
      for (int j=0;j<cont.size2();++j)
        PlusOp<T>()(cont(i,j),val(i,j),res(i,j));
  }
  else SpaceUnion(same,cont,val,res);
}

template<class T> void Catm(const genScalar<T> &cont, const genScalar<T> &val, genScalar<T> &res,const bool same[2])
{
   MinusOp<T>()(cont(),val(),res());
}

template<class T> void Catm(const genVector<T> &cont, const genVector<T> &val, genVector<T> &res,const bool same[2])
{
  if (same[0])
  {
    res.resize(cont.size(),false);
    for (int i=0;i<cont.size();++i)
      MinusOp<T>()(cont(i),val(i),res(i));
  }
  else SpaceUnionm(same,cont,val,res);
}

template<class T> void Catm(const genMatrix<T> &cont, const genMatrix<T> &val, genMatrix<T> &res,const bool same[2])
{
  if ((same[0])&&(same[1]))
  {
    res.resize(cont.size1(), cont.size2(),false);
    for (int i=0;i<cont.size1();++i)
      for (int j=0;j<cont.size2();++j)
        MinusOp<T>()(cont(i,j),val(i,j),res(i,j));
  }
  else SpaceUnionm(same,cont,val,res);
}


#endif // _GEN_MATRICIAL_TRAITS_H