// 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.666)


#ifndef _GEN_OPERATORS__H_
#define _GEN_OPERATORS__H_

#include "genTerm.h"
#include "genTermApply.h"
#include "genTermCompose.h"
#include "genTensorBuild.h"
#include "genTensorOps.h"
#include "savedGenTerm.h"


// sensible user interface
// important operators

template<class T,class GT>
std::shared_ptr<genTerm<T,GT::Nsp > > 
  Build(const std::shared_ptr<GT> &gt);

template<class T,class GT1,class GT2>
std::shared_ptr<genTerm<T,GT1::Nsp> > 
  Build(const std::shared_ptr<GT1> &gt1,const std::shared_ptr<GT2> &gt2);

template<class GT1,class GT2>
std::shared_ptr<genTerm<typename GT1::ValType,GT1::Nsp> >
  operator+(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2);

template<class GT>
std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
  operator*(const std::shared_ptr<GT> &gt, const typename GT::Scalar d);

template<class GT>
std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
  operator*(const typename GT::Scalar d, const std::shared_ptr<GT> &gt);

template<class GT>
std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
  Inverse(const std::shared_ptr<GT> &gt);

template<class GT1,class GT2>
std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::ContrProdType,GT1::Nsp> >
  Divide(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2);

template<class GT> 
std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> > 
  Transpose(const std::shared_ptr<GT> &gt, const int nn=0, const int mm=1);

template<class GT> 
std::shared_ptr<genTerm<typename GT::GradType,GT::Nsp> >
  Gradient(const std::shared_ptr<GT> &gt);

template<class GT>
std::shared_ptr<genTerm<typename GT::HessType,GT::Nsp> >
  Hessian(const std::shared_ptr<GT> &gt);

template<class GT>
std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
  Conj(const std::shared_ptr<GT> &gt);

template<class GT> 
std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
  Re(const std::shared_ptr<GT> &gt);

template<class GT> 
std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
  Im(const std::shared_ptr<GT> &gt);

template<class GT> 
std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
  Abs(const std::shared_ptr<GT> &gt);

template<class GT> 
std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
  Arg(const std::shared_ptr<GT> &gt);

  // blob of implementation details

template<class T,class GT> std::shared_ptr<genTerm<T,GT::Nsp > >
Build(const std::shared_ptr<GT> &gt)
{
  return Apply(gt,BuildFromOneOp<typename GT::ValType,T>());
}

template<class T,class GT1,class GT2> std::shared_ptr<genTerm<T,GT1::Nsp> >
Build(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Apply(gt1,gt2, BuildFromTwoOp<typename GT1::ValType,typename GT2::ValType,T>());
}


template<class GT> std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
operator*(const std::shared_ptr<GT> &gt, const typename GT::Scalar d)
{
  typename genTerm<genTensor0<typename GT::Scalar>,0>::Handle dterm(new ConstantField<genTensor0<typename GT::Scalar> >(d));
  return Compose<FullContrProdOp>(gt,dterm);
}


template<class GT> std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
operator*(const typename GT::Scalar d, const std::shared_ptr<GT> &gt)
{
  typename genTerm<genTensor0<typename GT::Scalar>,0>::Handle dterm(new ConstantField<genTensor0<typename GT::Scalar> >(d));
  return Compose<FullContrProdOp>(gt,dterm);
}

// new genTensor

// conflicts with default definition above
// need to specialize over genTerm<genTensor<> >
 
template<int order,class scalar, int N, int Nsp > std::shared_ptr<genTerm<genTensor<order,scalar,N>,Nsp> >
operator*(const std::shared_ptr<genTerm<genTensor<order,scalar,N>,Nsp> > &gt, const scalar d)
{
  typename genTerm<genTensor<0,scalar>,0>::Handle dterm(new ConstantField<genTensor<0,scalar> >(d));
  return Compose<FullContrProdOp>(gt,dterm);
}


template<int order,class scalar, int N, int Nsp > std::shared_ptr<genTerm<genTensor<order,scalar,N>,Nsp> >
operator*(const scalar d, const std::shared_ptr<genTerm<genTensor<order,scalar,N>,Nsp> > &gt)
{
  typename genTerm<genTensor<0,scalar>,0>::Handle dterm(new ConstantField<genTensor<0,scalar> >(d));
  return Compose<FullContrProdOp>(gt,dterm);
}


template<class GT> std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
Inverse (const std::shared_ptr<GT> &gt)
{
  return Apply<InverseOp>(gt);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::ContrProdType,GT1::Nsp> >
Divide (const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Apply<ContrProdOp>(gt1,Inverse(gt2));
}

template<class GT> std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> > 
Transpose(const std::shared_ptr<GT> &gt, const int nn, const int mm)
{
  return Apply(gt,TransposeOp<typename GT::ValType>(nn,mm));
}

template<class GT> std::shared_ptr<genTerm<typename GT::ValType,GT::Nsp> >
Conj(const std::shared_ptr<GT> &gt)
{
  return Apply<ConjOp>(gt);
}

template<class GT> std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
Re(const std::shared_ptr<GT> &gt)
{
  return Apply<ReOp>(gt);
}

template<class GT> std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
Im(const std::shared_ptr<GT> &gt)
{
  return Apply<ImOp>(gt);
}

template<class GT> std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
Abs(const std::shared_ptr<GT> &gt)
{
  return Apply<AbsOp>(gt);
}

template<class GT> std::shared_ptr<genTerm<typename TensorialTraits<typename GT::ValType>::RealValType,GT::Nsp> >
Arg(const std::shared_ptr<GT> &gt)
{
  return Apply<ArgOp>(gt);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::TensProdType,GT1::Nsp > >
TensorialProduct (const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Apply<TensProdOp>(gt1,gt2);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::FullContrProdType,GT1::Nsp > >
FullContractedProduct (const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Apply<FullContrProdOp>(gt1,gt2);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::TensProdType,GT1::Nsp+GT2::Nsp> >
TensorialProductCompose(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Compose<TensProdOp>(gt1,gt2);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::FullContrProdType,GT1::Nsp+GT2::Nsp> >
FullContractedProductCompose(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2 > &gt2)
{
  return Compose<FullContrProdOp>(gt1,gt2);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::ContrProdType,GT1::Nsp+GT2::Nsp> >
ContractedProductCompose(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Compose<ContrProdOp>(gt1,gt2);
}

template<class GT1,class GT2> std::shared_ptr<genTerm<typename TensorialTraitsBinary<typename GT1::ValType,typename GT2::ValType>::ContrProdType,GT1::Nsp+GT2::Nsp> >
NewContractedProductCompose(const std::shared_ptr<GT1> &gt1, const std::shared_ptr<GT2> &gt2)
{
  return Compose<NewContrProdOp>(gt1,gt2);
}

//--------------------------------------------------------------------------
// Non standard Operators
//--------------------------------------------------------------------------

#include "genPlusOperator.h"

template<class T, int n> class GradTerm : public genTerm< typename TensorialTraits<T>::GradType, n>
{
public :
  typedef typename TensorialTraits<T>::Scalar Scalar;
  typedef typename TensorialTraits<T>::GradType ValType;
  typedef typename ContainerTraits<ValType,n>::Container ContainerValType;
  static const int Nsp=n;
protected :
  typename diffTerm<T,n>::ConstHandle space1;
  
public :
  GradTerm(const typename diffTerm<T,n>::ConstHandle &space1_) : space1(space1_) {}
  virtual int getNumKeys(MElement* ele,int k=0) const {return space1->getNumKeys(ele,k);}
  virtual void getKeys(MElement* ele, std::vector<Dof> &keys, int k=0) const { space1->getKeys(ele,keys,k);}
  virtual int getIncidentSpaceTag(int k=0) const { return space1->getIncidentSpaceTag(k);}
  virtual void get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const;
//  virtual void get(MElement* ele, int npts, IntPt* GP, ContainerValType &vals) const {LinearTermBase<typename TensorialTraits<T1>::ValType>::get(ele,npts,GP,vals);}
  virtual GradTerm<T,n>* clone () const { return new GradTerm<T,n>(space1);}
  virtual ~GradTerm() {}
};

template<class GT> std::shared_ptr<genTerm<typename GT::GradType,GT::Nsp> >
Gradient(const std::shared_ptr<GT> &gt)
{
  return std::shared_ptr<genTerm<typename GT::GradType,GT::Nsp> >(new GradTerm<typename GT::ValType,GT::Nsp>(gt));
}



template<class T, class T2, int n> class GradTerm2 : public genTerm< T2, n>
{
public :
  typedef typename TensorialTraits<T>::Scalar Scalar;
  typedef T2 ValType;
  typedef typename ContainerTraits<ValType,n>::Container ContainerValType;
  static const int Nsp=n;
protected :
  typename diffTerm<T,n>::ConstHandle space1;
  
public :
  GradTerm2(const typename diffTerm<T,n>::ConstHandle &space1_) : space1(space1_) {}
  virtual int getNumKeys(MElement* ele,int k=0) const {return space1->getNumKeys(ele,k);}
  virtual void getKeys(MElement* ele, std::vector<Dof> &keys, int k=0) const { space1->getKeys(ele,keys,k);}
  virtual int getIncidentSpaceTag(int k=0) const { return space1->getIncidentSpaceTag(k);}
  virtual void get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const;
//  virtual void get(MElement* ele, int npts, IntPt* GP, ContainerValType &vals) const {LinearTermBase<typename TensorialTraits<T1>::ValType>::get(ele,npts,GP,vals);}
  virtual GradTerm2<T,T2,n>* clone () const { return new GradTerm2<T,T2,n>(space1);}
  virtual ~GradTerm2() {}
};

template<class GT,class T2> std::shared_ptr<genTerm<T2,GT::Nsp> >
Gradient(const std::shared_ptr<GT> &gt)
{
  return std::shared_ptr<genTerm<T2,GT::Nsp> >(new GradTerm2<typename GT::ValType,T2,GT::Nsp>(gt));
}



template<class T, int n> class HessTerm : public genTerm< typename TensorialTraits<T>::HessType, n>
{
public :
  typedef typename TensorialTraits<T>::Scalar Scalar;
  typedef typename TensorialTraits<T>::HessType ValType;
  typedef typename ContainerTraits<ValType,n>::Container ContainerValType;
  static const int Nsp=n;
protected :
  typename diffTerm<T,n>::ConstHandle space1;
  
public :
  HessTerm(const typename diffTerm<T,n>::ConstHandle &space1_) : space1(space1_) {}
  virtual int getNumKeys(MElement* ele,int k=0) const {return space1->getNumKeys(ele,k);}
  virtual void getKeys(MElement* ele, std::vector<Dof> &keys, int k=0) const { space1->getKeys(ele,keys,k);}
  virtual int getIncidentSpaceTag(int k=0) const { return space1->getIncidentSpaceTag(k);}
  virtual void get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const;
//  virtual void get(MElement* ele, int npts, IntPt* GP, ContainerValType &vals) const {LinearTermBase<typename TensorialTraits<T1>::ValType>::get(ele,npts,GP,vals);}
  virtual HessTerm<T,n>* clone () const { return new HessTerm<T,n>(space1);}
  virtual ~HessTerm() {}
};

template<class GT> std::shared_ptr<genTerm<typename GT::HessType,GT::Nsp> >
Hessian(const std::shared_ptr<GT> &gt)
{
  return std::shared_ptr<genTerm<typename GT::HessType,GT::Nsp> >(new HessTerm<typename GT::ValType,GT::Nsp>(gt));
}



//--------------------------------------------------------------------------
// Get
//--------------------------------------------------------------------------

template<class T,int n> void GradTerm<T,n>::get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const
{
  space1->getgradf(ele,npts,GP,vvals);
}

template<class T,class T2,int n> void GradTerm2<T,T2,n>::get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const
{
  space1->getgradf(ele,npts,GP,vvals);
}

template<class T,int n> void HessTerm<T,n>::get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const
{
  space1->gethessf(ele,npts,GP,vvals);
}


#endif // _GEN_OPERATORS_H_