gsMaterialMatrixIntegrate.h

Go to the documentation of this file.

#pragma once

#include <gsCore/gsFunction.h>
#include <gsKLShell/src/gsMaterialMatrixContainer.h>

namespace gismo
{

template <class T, enum MaterialOutput out> class gsMaterialMatrixIntegrateSingle;

template <class T, enum MaterialOutput out>
class gsMaterialMatrixIntegrate : public gsFunction<T>
{
public:
    gsMaterialMatrixIntegrate(  const gsMaterialMatrixContainer<T> & materialMatrices,
                            const gsFunctionSet<T> * deformed)
    :
    m_materialMatrices(materialMatrices),
    m_deformed(deformed)
    {
        for (index_t p = 0; p!=m_materialMatrices.size(); p++)
        {
            GISMO_ASSERT(materialMatrices.piece(p)!=nullptr,"Material matrix "<<p<<" is incomplete!");
            GISMO_ASSERT(materialMatrices.piece(p)!=NULL,"Material matrix "<<p<<" is incomplete!");
            GISMO_ASSERT(materialMatrices.piece(p)->initialized(),"Material matrix "<<p<<" is incomplete!");
        }

        this->_makePieces(deformed);
    }

    gsMaterialMatrixIntegrate(  gsMaterialMatrixBase<T> * materialMatrix,
                            const gsFunctionSet<T> * deformed)
    :
    m_materialMatrices(deformed->nPieces()),
    m_deformed(deformed)
    {
        GISMO_ASSERT(materialMatrix->initialized(),"Material matrix is incomplete!");
        for (index_t p = 0; p!=deformed->nPieces(); ++p)
            m_materialMatrices.set(p,materialMatrix);
        this->_makePieces(deformed);
    }

    gsMaterialMatrixIntegrate(  const gsMaterialMatrixContainer<T> & materialMatrices,
                            const gsFunctionSet<T> * undeformed,
                            const gsFunctionSet<T> * deformed)
    :
    m_materialMatrices(materialMatrices),
    m_deformed(deformed)
    {
        this->_makePieces(undeformed,deformed);
    }

    gsMaterialMatrixIntegrate(  gsMaterialMatrixBase<T> * materialMatrix,
                            const gsFunctionSet<T> * undeformed,
                            const gsFunctionSet<T> * deformed)
    :
    m_materialMatrices(deformed->nPieces()),
    m_deformed(deformed)
    {
        for (index_t p = 0; p!=deformed->nPieces(); ++p)
        {
            m_materialMatrices.set(p,materialMatrix);
            this->_makePieces(undeformed,deformed);
        }
    }

    ~gsMaterialMatrixIntegrate()
    {
        freeAll(m_pieces);
    }

    short_t domainDim() const {return 2;}

    short_t targetDim() const { return this->piece(0).targetDim(); }

    const gsFunction<T> & piece(const index_t p) const
    {
        return *m_pieces[p];
    }

    void eval_into(const gsMatrix<T>& u, gsMatrix<T>& result) const
    { GISMO_NO_IMPLEMENTATION; }

protected:
    void _makePieces(const gsFunctionSet<T> * deformed)
    {
        m_pieces.resize(deformed->nPieces());
        for (size_t p = 0; p!=m_pieces.size(); ++p)
            m_pieces.at(p) = new gsMaterialMatrixIntegrateSingle<T,out>(p,m_materialMatrices.piece(p),deformed);
    }

    void _makePieces(const gsFunctionSet<T> * undeformed, const gsFunctionSet<T> * deformed)
    {
        m_pieces.resize(m_deformed->nPieces());
        for (size_t p = 0; p!=m_pieces.size(); ++p)
            m_pieces.at(p) = new gsMaterialMatrixIntegrateSingle<T,out>(p,m_materialMatrices.piece(p),undeformed,deformed);
    }

protected:
    gsMaterialMatrixContainer<T> m_materialMatrices;
    const gsFunctionSet<T> * m_deformed;
    mutable std::vector<gsMaterialMatrixIntegrateSingle<T,out> *> m_pieces;
};

template <class T, enum MaterialOutput out>
class gsMaterialMatrixIntegrateSingle : public gsFunction<T>
{
public:

    gsMaterialMatrixIntegrateSingle(  index_t patch,
                                gsMaterialMatrixBase<T> * materialMatrix,
                                const gsFunctionSet<T> * deformed);

    gsMaterialMatrixIntegrateSingle(  index_t patch,
                                gsMaterialMatrixBase<T> * materialMatrix,
                                const gsFunctionSet<T> * undeformed,
                                const gsFunctionSet<T> * deformed);

    gsMaterialMatrixIntegrateSingle(){};

    short_t domainDim() const {return 2;}

    short_t targetDim() const { return targetDim_impl<out>(); }


private:
    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::Density   , short_t>::type targetDim_impl() const { return 1; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorN ||
                            _out==MaterialOutput::CauchyVectorN ||
                            _out==MaterialOutput::VectorM ||
                            _out==MaterialOutput::CauchyVectorM   , short_t>::type targetDim_impl() const { return 3; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixA ||
                            _out==MaterialOutput::MatrixB ||
                            _out==MaterialOutput::MatrixC ||
                            _out==MaterialOutput::MatrixD   , short_t>::type targetDim_impl() const { return 9; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::PStressN ||
                            _out==MaterialOutput::PStressM ||
                            _out==MaterialOutput::PStrainN ||
                            _out==MaterialOutput::PStrainM   , short_t>::type targetDim_impl() const { return 2; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::Stretch   , short_t>::type targetDim_impl() const { return 3; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::StretchDir, short_t>::type targetDim_impl() const { return 9; };

protected:
    void setPatch(index_t p) {m_pIndex = p; }

public:
    void eval_into(const gsMatrix<T>& u, gsMatrix<T>& result) const;

private:
    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::Density   , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorN ||
                            _out==MaterialOutput::CauchyVectorN   , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorM ||
                            _out==MaterialOutput::CauchyVectorM   , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixA ||
                            _out==MaterialOutput::MatrixB ||
                            _out==MaterialOutput::MatrixC ||
                            _out==MaterialOutput::MatrixD   , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::PStressN ||
                            _out==MaterialOutput::PStressM  , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::PStrainN ||
                            _out==MaterialOutput::PStrainM  , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::Stretch   , void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::StretchDir, void>::type eval_into_impl(const gsMatrix<T>& u, gsMatrix<T>& result) const;

protected:

    T getMoment() const { return getMoment_impl<out>(); }

private:
    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorN ||
                            _out==MaterialOutput::PStressN||
                            _out==MaterialOutput::CauchyVectorN, T>::type getMoment_impl() const { return 0; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorM ||
                            _out==MaterialOutput::PStressM||
                            _out==MaterialOutput::CauchyVectorM, T>::type getMoment_impl() const { return 1; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixA, T>::type getMoment_impl() const { return 0; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixB, T>::type getMoment_impl() const { return 1; };

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixC, T>::type getMoment_impl() const { return 1; }; // must be 1

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixD, T>::type getMoment_impl() const { return 2; }; // must be 2

    template<enum MaterialOutput _out>
    typename std::enable_if<!(_out==MaterialOutput::VectorN         || _out==MaterialOutput::VectorM        ||
                              _out==MaterialOutput::CauchyVectorN   || _out==MaterialOutput::CauchyVectorM  ||
                              _out==MaterialOutput::MatrixA         || _out==MaterialOutput::MatrixB        ||
                              _out==MaterialOutput::MatrixC         || _out==MaterialOutput::MatrixD        ||
                              _out==MaterialOutput::PStressN        || _out==MaterialOutput::PStressM  )
                                                 , index_t>::type getMoment_impl() const { GISMO_NO_IMPLEMENTATION};

protected:

    void integrateZ_into(const gsMatrix<T> & u, const index_t moment, gsMatrix<T> & result) const;

    void multiplyLinZ_into(const gsMatrix<T> & u, const index_t moment, gsMatrix<T> & result) const;

    void multiplyZ_into(const gsMatrix<T> & u, index_t moment, gsMatrix<T> & result) const;

protected:

    gsMatrix<T> _eval3D(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;


    gsMatrix<T> _eval  (const gsMatrix<T>& u) const;

private:
    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::VectorN ||
                            _out==MaterialOutput::VectorM         , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::CauchyVectorN ||
                            _out==MaterialOutput::CauchyVectorM   , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;


    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::MatrixA ||
                            _out==MaterialOutput::MatrixB ||
                            _out==MaterialOutput::MatrixC ||
                            _out==MaterialOutput::MatrixD   , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<_out==MaterialOutput::PStress  ||
                            _out==MaterialOutput::PStressN ||
                            _out==MaterialOutput::PStressM   , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const;

    template<enum MaterialOutput _out>
    typename std::enable_if<!(_out==MaterialOutput::VectorN         || _out==MaterialOutput::VectorM        ||
                              _out==MaterialOutput::CauchyVectorN   || _out==MaterialOutput::CauchyVectorM  ||
                              _out==MaterialOutput::MatrixA         || _out==MaterialOutput::MatrixB        ||
                              _out==MaterialOutput::MatrixC         || _out==MaterialOutput::MatrixD        ||
                              _out==MaterialOutput::PStressN        || _out==MaterialOutput::PStressM  )
                                                            , gsMatrix<T>>::type eval3D_impl(const gsMatrix<T>& u, const gsMatrix<T>& Z) const
    { GISMO_NO_IMPLEMENTATION};

protected:
    index_t m_pIndex;
    gsMaterialMatrixBase<T> * m_materialMat;
    gsMatrix<T> m_z;
};

} // namespace


#ifndef GISMO_BUILD_LIB
#include GISMO_HPP_HEADER(gsMaterialMatrixIntegrate.hpp)
#endif