gsPreCICEFunction.h

Go to the documentation of this file.

#pragma once

#include <gsCore/gsLinearAlgebra.h>
#include <gsCore/gsGeometry.h>
#include <gsUtils/gsCombinatorics.h>
#include <gsPreCICE/gsPreCICE.h>

namespace gismo
{

template <class T>
class gsPreCICEFunction : public gsFunction<T>
{
public:
    typedef gsGeometry<T> Base;

    typedef memory::shared_ptr< gsPreCICEFunction > Ptr;

    typedef memory::unique_ptr< gsPreCICEFunction > uPtr;

    gsPreCICEFunction() { }

    gsPreCICEFunction(        gsPreCICE<T> *    interface,
                        const index_t &         meshID,
                        const index_t &         dataID,
                        const gsMultiPatch<T> & patches,
                        const bool parametric = false)
    :
    m_interface(interface),
    m_meshID(meshID),
    m_dataID(dataID),
    m_patches(patches),
    m_parametric(parametric),
    m_patchID(0),
    m_domainDim(m_patches.domainDim())
    {
    }

    static uPtr make(   const gsPreCICE<T> *    interface,
                        const index_t &         meshID,
                        const index_t &         dataID,
                        const gsMultiPatch<T> & patches,
                        const bool parametric = false)
    { return uPtr(new gsPreCICEFunction(interface, meshID, dataID, patches, parametric)); }

    GISMO_CLONE_FUNCTION(gsPreCICEFunction)

    
    const gsPreCICEFunction<T> & piece(const index_t) const
    {
        return *this;
    }

    virtual short_t domainDim() const
    { return m_domainDim; }

    virtual short_t targetDim() const
    { return 1; }

    virtual void eval_into(const gsMatrix<T>& u, gsMatrix<T>& result) const
    {
        GISMO_ASSERT(u.rows() == m_domainDim, "Wrong domain dimension "<< u.rows()
                                              << ", expected "<< m_domainDim);

        // Does not work
        // if (m_parametric)
        //     m_interface->readBlockScalarData(m_meshID,m_dataID,m_patches.patch(m_patchID).eval(u),result);
        // if (m_parametric)
        //     m_interface->readBlockScalarData(m_meshID,m_dataID,u,result);

        gsMatrix<T> coords(m_patches.targetDim(),u.cols());
        if (m_parametric)
            m_patches.patch(m_patchID).eval_into(u,coords);
        else
            coords = u;

        m_interface->readBlockScalarData(m_meshID,m_dataID,coords,result);
    }

    virtual void deriv_into(const gsMatrix<T>& u, gsMatrix<T>& result) const
    {
        // This would be nice to have with higher-order (IGA) coupling of precice
        GISMO_NO_IMPLEMENTATION;
    }

    virtual void deriv2_into(const gsMatrix<T>& u, gsMatrix<T>& result) const
    {
        // This would be nice to have with higher-order (IGA) coupling of precice
        GISMO_NO_IMPLEMENTATION;
    }

    void evalAllDers_into(const gsMatrix<T> & u, int n,
                          std::vector<gsMatrix<T> > & result) const
    {
        // This would be nice to have with higher-order (IGA) coupling of precice
        gsMatrix<T> tmp;
        this->eval_into(u,tmp);
        result.push_back(tmp);
    }

    virtual std::ostream &print(std::ostream &os) const
    {
        os << "gsPreCICEFunction, defined using gsPreCICE (precice::SolverInterface):\n";
        os << m_interface;
        return os;
    }


private:

    gsPreCICE<T> * m_interface;
    index_t m_meshID, m_dataID;
    gsMultiPatch<T> m_patches;
    bool m_parametric;
    index_t m_patchID;
    index_t m_domainDim;

};

}