// 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: Eric Bechet (rev.6987 in gmsh)


#ifndef _GENSOLVERFIELD_H_
#define _GENSOLVERFIELD_H_

#include <vector>
#include "MElement.h"
#include "genTerm.h"
#include "dofManager.h"
#include "genFSpace.h"


// this is a field that is the result of a FE computation
// values of dofs are known and stored in a dofManager.

template<class T>
class genSField : public diffTerm<T,0>
{
public:
  typedef typename TensorialTraits<T>::Scalar Scalar;
  typedef typename TensorialTraits<T>::ValType ValType;
  typedef typename TensorialTraits<T>::GradType GradType;
  typedef typename TensorialTraits<T>::HessType HessType;
  typedef typename ContainerTraits<ValType,0>::Container ContainerValType;
  typedef typename ContainerTraits<GradType,0>::Container ContainerGradType;
  typedef typename ContainerTraits<HessType,0>::Container ContainerHessType;
  static const int Nsp=0;
private:
  dofManager<Scalar>* dm;
  typename diffTerm<T,1>::ConstHandle fs;

public:
  genSField(dofManager<Scalar>* dm_, const typename diffTerm<T,1>::ConstHandle &fs_) : dm(dm_), fs(fs_) {}
  virtual int getNumKeys(MElement* ele,int k=0) const  {return 0;}
  virtual void getKeys(MElement* ele, std::vector<Dof> &keys, int k=0) const {}
  virtual int getIncidentSpaceTag(int k=0) const { return 0;}

  virtual void get(MElement* ele, int npts, IntPt* GP, std::vector<ContainerValType> &vvals) const
  {
    std::vector<Dof> D;
    std::vector<Scalar> DMVals;
    std::vector< typename genFSpace<T>::ContainerValType > vFSvals(npts);

    fs->getKeys(ele, D);
    dm->getDofValue(D, DMVals);
    fs->get(ele, npts, GP, vFSvals);

    int nbdofs = D.size();
    int curpos = 0;
    for (unsigned int i = 0; i < npts; ++i)
    {
      vvals[i]=ContainerValType();
      for (unsigned int j = 0; j < nbdofs; ++j)
        vvals[i]() += vFSvals[i](j) * genTensor0<Scalar>(DMVals[j])();
    }
  }
  virtual void get(MElement* ele, int npts, IntPt* GP, ContainerValType &vals) const
  {
    diffTerm<T,0>::get(ele,npts,GP,vals);
  }
  virtual void getgradf(MElement* ele, int npts, IntPt* GP, std::vector< ContainerGradType > &vgrads) const
  {
    std::vector<Dof> D;
    std::vector<Scalar> DMVals;
    std::vector< typename genFSpace<T>::ContainerGradType> vSFGrads(npts);

    fs->getKeys(ele, D);
    dm->getDofValue(D, DMVals);
    fs->getgradf(ele, npts, GP, vSFGrads);

    int nbdofs = D.size();
    int curpos = 0;

    for (unsigned int i = 0; i < npts; ++i)
    {
      vgrads[i]=ContainerGradType();
      for (unsigned int j = 0; j < nbdofs; ++j)
        vgrads[i]() += vSFGrads[i](j) * genTensor0<Scalar>(DMVals[j])();
    }
  }
  virtual void getgradf(MElement* ele, int npts, IntPt* GP, ContainerGradType &grads) const
  {
    diffTerm<T,0>::getgradf(ele,npts,GP,grads);
  }
  virtual void getgradfuvw(MElement* ele, int npts, IntPt* GP, std::vector< ContainerGradType > &vgrads) const
  {
    std::vector<Dof> D;
    std::vector<Scalar> DMVals;
    std::vector< typename genFSpace<T>::ContainerGradType> vSFGrads(npts);

    fs->getKeys(ele, D);
    dm->getDofValue(D, DMVals);
    fs->getgradfuvw(ele, npts, GP, vSFGrads);

    int nbdofs = D.size();
    int curpos = 0;

    for (unsigned int i = 0; i < npts; ++i)
    {
      vgrads[i]=ContainerGradType();
      for (unsigned int j = 0; j < nbdofs; ++j)
        vgrads[i]() += vSFGrads[i](j) * genTensor0<Scalar>(DMVals[j])();
    }
  }
  virtual void gethessf(MElement* ele, int npts, IntPt* GP, std::vector< ContainerHessType > &vhesss) const
  {
    std::vector<Dof> D;
    std::vector<Scalar> DMVals;
    std::vector< typename genFSpace<T>::ContainerHessType> vSFHesss(npts);

    fs->getKeys(ele, D);
    dm->getDofValue(D, DMVals);
    fs->gethessf(ele, npts, GP, vSFHesss);

    int nbdofs = D.size();
    int curpos = 0;

    for (unsigned int i = 0; i < npts; ++i)
    {
      vhesss[i]=ContainerHessType();
      for (unsigned int j = 0; j < nbdofs; ++j)
        vhesss[i]() += vSFHesss[i](j) * genTensor0<Scalar>(DMVals[j])();
    }
  }
  virtual void gethessf(MElement* ele, int npts, IntPt* GP, ContainerHessType &hesss) const
  {
    diffTerm<T,0>::gethessf(ele,npts,GP,hesss);
  }

  virtual genSField<T>* clone () const { return new genSField<T>(dm,fs);}
};

#endif //_GENSOLVERFIELD_H_