gsPatchRule.h

Go to the documentation of this file.

 
#pragma once
 
#include <gsAssembler/gsQuadRule.h>
 
namespace gismo
{
 
template<class T>
class gsPatchRule GISMO_FINAL : public gsQuadRule<T>
{
public:
 
    typedef memory::unique_ptr<gsPatchRule> uPtr;
 
    gsPatchRule()
    :
    m_basis(nullptr),
    m_deg(0),
    m_reg(0),
    m_over(false),
    m_fixDir(-1)
    {};
 
    gsPatchRule(const gsBasis<T> & basis,
                const index_t degree,
                const index_t regularity,
                const bool overintegrate,
                const short_t fixDir = -1);
 
 
    static uPtr make(   const gsBasis<T> & basis,
                        const index_t degree,
                        const index_t regularity,
                        const bool overintegrate,
                        const short_t fixDir = -1)
    { return uPtr( new gsPatchRule(basis,degree,regularity,overintegrate,fixDir) ); }
 
 
    ~gsPatchRule() { };
 
public:
 
    void mapTo( const gsVector<T>& lower, const gsVector<T>& upper,
                       gsMatrix<T> & nodes, gsVector<T> & weights ) const;
 
    void mapTo(  T startVal, T endVal,
                       gsMatrix<T> & nodes, gsVector<T> & weights ) const;
 
    void mapToAll( const std::vector<T> &,
                   gsMatrix<T> &, gsVector<T> &) const
    { GISMO_NO_IMPLEMENTATION }
 
    index_t dim() const { return m_basis->dim(); }
 
    gsVector<T> nodes(const index_t d=0) const { return m_nodes.at(d); }
    gsVector<T> weights(const index_t d=0) const { return m_weights.at(d); }
 
protected:
 
    index_t _numQuads(const gsKnotVector<T> knots) const
    {
        index_t ndof = (m_deg + 1)*knots.numElements() - (m_reg+1)*(knots.numElements()-1);
        ndof = static_cast<index_t>(std::ceil(ndof/2.0));
        return ndof;
    }
 
    bool _isSymmetric(const gsKnotVector<T> knots) const
    {
        gsKnotVector<T> copy = knots;
        copy.reverse();
        std::vector<T> result;
        knots.symDifference(copy,result);
        return result.size()==0 ? true : false;
    }
 
 
    gsKnotVector<T> _init(const gsBSplineBasis<T> * Bbasis) const;
 
    std::pair<gsMatrix<T>,gsVector<T> > _integrate(const gsKnotVector<T> & knots ) const;
 
    std::pair<gsVector<T>,gsVector<T> > _compute(const gsKnotVector<T> & knots,
                                                const gsMatrix<T> & greville,
                                                const gsVector<T> & integrals,
                                                const T tol = 1e-10) const;
 
private:
    const gsBasis<T> * m_basis;
    const index_t m_deg,m_reg;
    const bool m_over;
    const short_t m_fixDir;
 
    std::vector<gsVector<T> > m_nodes;
    std::vector<gsVector<T> > m_weights;
 
    gsVector<T> m_end;
 
    mutable typename gsSparseSolver<T>::QR m_solver;
 
    size_t m_dim;
 
    std::vector<std::map<T,T> > m_maps;
}; // class gsPatchRule
 
} // namespace gismo
 
 
#ifndef GISMO_BUILD_LIB
#include GISMO_HPP_HEADER(gsPatchRule.hpp)
#endif