discreteFunctionShellFE.h

Go to the documentation of this file.

#ifndef __DISCRETEFUNCTIONSHELLFE_H
#define __DISCRETEFUNCTIONSHELLFE_H

#include <EigenDataContainer.h>
#include <initializer_list>

namespace shellFE {

template <typename ConfiguratorType, typename VectorType, typename NLTYPE >
struct DiscreteFunctionLookup {

  void getLocalDof ( const ConfiguratorType & conf, const VectorType &dofs, const typename ConfiguratorType::ElementType &El,
                     const int bfNum,  const NLTYPE &bfValue, NLTYPE &aux ) {
      aux = bfValue;
      aux *= dofs[ conf.localToGlobal ( El, bfNum ) ];
  }

};

template <typename ConfiguratorType, typename VectorType = typename ConfiguratorType::VectorType >
class DiscreteFunctionDefaultShellFE {
public:
  typedef typename ConfiguratorType::RealType RealType;
  typedef typename ConfiguratorType::DomVecType  DomVecType;
  typedef typename ConfiguratorType::TangentVecType TangentVecType;
  typedef typename ConfiguratorType::Matrix22  Matrix22;
  typedef typename ConfiguratorType::ElementType ElementType;
  typedef typename ConfiguratorType::BaseFuncSetType BaseFuncSetType;

  const ConfiguratorType & _conf;
  const VectorType & _dofs;

  DiscreteFunctionDefaultShellFE ( const ConfiguratorType & config, const VectorType & Dofs )
    : _conf ( config ),
      _dofs ( Dofs ) { }

  RealType evaluate ( const ElementType &El, const DomVecType& RefCoord ) const {
    const BaseFuncSetType &bfs = _conf.getBaseFunctionSet ( El );
    RealType w = 0.;
    RealType v, aux;
    for ( int b = 0; b < static_cast<int> ( _conf.getNumLocalDofs ( El ) ); ++b ) {
      v = bfs.evaluate ( b, RefCoord );
      DiscreteFunctionLookup<ConfiguratorType, VectorType, RealType>().getLocalDof( _conf, _dofs, El, b, v, aux );
      w += aux;
    }
    return w;
  }

  RealType evaluateAtQuadPoint ( const ElementType &El, int QuadPoint ) const {
    const BaseFuncSetType &bfs = _conf.getBaseFunctionSet ( El );
    RealType w = 0.;
    RealType v, aux;
    for ( int b = 0; b < static_cast<int> ( _conf.getNumLocalDofs ( El ) ); ++b ) {
      v = bfs.evaluate ( b, QuadPoint );
      DiscreteFunctionLookup<ConfiguratorType, VectorType, RealType>().getLocalDof ( _conf, _dofs, El, b, v, aux );
      w += aux;
    }
    return w;
  }

  void evaluateGradient ( const ElementType &El, const DomVecType& RefCoord, DomVecType& Grad ) const {
    const BaseFuncSetType &bfs = _conf.getBaseFunctionSet ( El );
    DomVecType v, aux;
    Grad.setZero();
    for ( int b = 0; b < static_cast<int> ( _conf.getNumLocalDofs ( El ) ); ++b ) {
      bfs.evaluateGradient ( b, RefCoord, v );
      DiscreteFunctionLookup<ConfiguratorType, VectorType, DomVecType>().getLocalDof ( _conf, _dofs, El, b, v, aux );
      Grad += aux;
    }
  }

  void evaluateGradientAtQuadPoint ( const ElementType &El, int QuadPoint, DomVecType& Grad ) const {
    const BaseFuncSetType &bfs = _conf.getBaseFunctionSet ( El );
    DomVecType v, aux;
    Grad.setZero();
    for ( int b = 0; b < static_cast<int> ( _conf.getNumLocalDofs ( El ) ); ++b ) {
      v = bfs.evaluateGradient ( b, QuadPoint );
      DiscreteFunctionLookup<ConfiguratorType, VectorType, DomVecType >().getLocalDof ( _conf, _dofs, El, b, v, aux );
      Grad += aux;
    }
  }

  const VectorType& getDofs (  ) const {return _dofs;}

};

template <typename ConfiguratorType>
class DiscreteVectorFunctionDefaultShellFE {
public:
  typedef typename ConfiguratorType::RealType RealType;
  typedef typename ConfiguratorType::DomVecType  DomVecType;
  typedef typename ConfiguratorType::TangentVecType TangentVecType;
  typedef typename ConfiguratorType::VectorType VectorType;
  typedef typename ConfiguratorType::Matrix22  Matrix22;
  typedef typename ConfiguratorType::Matrix32  Matrix32;
  typedef typename ConfiguratorType::Matrix33  Matrix33;

  typedef typename ConfiguratorType::ElementType ElementType;
  typedef DiscreteFunctionDefaultShellFE<ConfiguratorType, Eigen::Ref<const VectorType> > DiscFuncType;

  const int _numComponents;
  mutable std::vector<Eigen::Ref<const VectorType> > _refs;
  mutable std::vector<DiscFuncType> _discrFuncs;

public:

  DiscreteVectorFunctionDefaultShellFE ( const ConfiguratorType &Configurator,  const VectorType &Dofs, const int numComponents ) :
      _numComponents ( numComponents )
  {
     const int numGlobalDofs = Configurator.getNumGlobalDofs();
     _refs.reserve( _numComponents );
     _discrFuncs.reserve( _numComponents );

     _refs.push_back ( Dofs.segment( 0, numGlobalDofs) );
     _refs.push_back ( Dofs.segment( numGlobalDofs, numGlobalDofs) );
     _refs.push_back ( Dofs.segment( 2 * numGlobalDofs, numGlobalDofs) );

     _discrFuncs.push_back ( DiscFuncType ( Configurator, _refs[0] ) );
     _discrFuncs.push_back ( DiscFuncType ( Configurator, _refs[1] ) );
     _discrFuncs.push_back ( DiscFuncType ( Configurator, _refs[2] ) );
  }

  void evaluate( const ElementType &El, const DomVecType& RefCoord, TangentVecType &Value ) const {
    for ( int c = 0; c < 3; ++c )
      Value[c] = _discrFuncs[c].evaluate ( El, RefCoord );
  }

  void evaluateAtQuadPoint ( const ElementType &El, int QuadPoint, TangentVecType &Value ) const {
    for ( int c = 0; c < 3; ++c )
      Value[c] = _discrFuncs[c].evaluateAtQuadPoint ( El, QuadPoint );
  }

  void evaluateGradient ( const ElementType &El, const DomVecType& RefCoord, Matrix32 &Dx ) const {
    DomVecType v;
    for ( int c = 0; c < 3; c++ ) {
      _discrFuncs[c].evaluateGradient ( El, RefCoord, v );
      Dx( c, 0 ) = v[0];
      Dx( c, 1 ) = v[1];
    }
  }

  void evaluateGradientAtQuadPoint ( const ElementType &El, int QuadPoint, Matrix32 &Dx ) const {
    DomVecType v;
    for ( int c = 0; c < 3; ++c ) {
      _discrFuncs[c].evaluateGradientAtQuadPoint ( El, QuadPoint, v );
      Dx( c, 0 ) = v[0];
      Dx( c, 1 ) = v[1];
    }
  }

  void evaluateFirstFundamentalForm ( const ElementType &El, const DomVecType& RefCoord, Matrix22 &g ) const {
    Matrix32 Dx;
    evaluateGradient( El, RefCoord, Dx );
    g = Dx.transpose() * Dx;
  }


  void evaluateFirstFundamentalFormAtQuadPoint ( const ElementType &El, int QuadPoint, Matrix22 &g ) const {
    Matrix32 Dx;
    evaluateGradientAtQuadPoint( El, QuadPoint, Dx );
    g = Dx.transpose() * Dx;
  }

  void evaluateFirstFundamentalForm ( const Matrix32 &Dx, Matrix22 &g ) const {
    g = Dx.transpose() * Dx;
  }

  // given chart X -> compute n \circ X  = d_1 X \cross d_2 X / |.|
  void evaluateNormalAtQuadPointOnShell ( const ElementType &El, int QuadPoint, TangentVecType &normal ) const {
    Matrix32 Dx;
    evaluateGradientAtQuadPoint ( El, QuadPoint, Dx );
    TangentVecType unNormalizedNormal ( (Dx.col(0)).cross(Dx.col(1)) );
    normal = unNormalizedNormal.normalized();
  }

  void evaluateNormal ( const Matrix32 &Dx, TangentVecType &normal ) const {
    TangentVecType unNormalizedNormal ( (Dx.col(0)).cross(Dx.col(1)) );
    normal = unNormalizedNormal.normalized();
  }


  const DiscFuncType& operator[] ( int i ) const { return _discrFuncs[i];}
  DiscFuncType& operator[] ( int i ) { return _discrFuncs[i];}

};

} // end namespace

#endif