NetElement3N.h

 
 
#ifndef C_NET_ELEMENT_3N
#define C_NET_ELEMENT_3N
 
#include <math.h>
#include "NetElementDef.h"
#include "NetElement3NForces.h"
#include <fhsim_environment/hydrodynamics/NetSolidity.h>
#include <fhsim/PrintDuringExec.h>
#include <sfh/ad/Dual.h>
 
#ifdef FH_VISUALIZATION
    //#define USE_DYNAMIC_LINES
    #ifdef USE_DYNAMIC_LINES
        #include "sfh/ogre/CDynamicLines.h"
    #else
        #include "sfh/ogre/C3DLine.h"
    #endif
#endif
 
#define RENDER_OUTLINE
#define RENDER_MESH
//#define ATTACH_SPHERE
//#define RENDER_MATERIAL
//#define FILTER_FORCES
//#define CNETELEMENT3N_TRANSPARENT
 
class CDynamicLines; // forward declaration
 
struct NetPanelOptions
{
    bool   twineStructure = true;
    double dampingRatio   = 1.0; 
};
 
class NetElement3N
{
public:
    NetElement3N(
        const std::string       netName,
        const unsigned long     iID,
        const NetElementDef*    netPanel,
        const double            dRho_water,
        const double            knotDiameter,
        const net_solidity::PanelSolidity& solidity,
        const NetPanelOptions&  options
        );
 
    ~NetElement3N(void);
 
    // Owns raw visualization buffers (freed in the destructor), so it must not be copied.
    NetElement3N(const NetElement3N&) = delete;
    NetElement3N& operator=(const NetElement3N&) = delete;
 
    template<class T>
    void AddNodeForces(const hydrodynamics::PanelLoadLaw& law, const hydrodynamics::Fluid& fluid,
        const T* const posA_ned, const T* const posB_ned, const T* const posC_ned,
        const T* const velA_ned, const T* const velB_ned, const T* const velC_ned,
        const T* const waterVel_ned,
        T* const nodeAforce_ned, T* const nodeBforce_ned, T* const nodeCforce_ned,
        double addedLinearDrag = 0.0, const hydrodynamics::PanelClampLog* clampLog = nullptr) const
    {
        const net_element_forces::PanelFrame<T> frame =
            net_element_forces::ComputePanelFrame(m_panelForceParams, posA_ned, posB_ned, posC_ned);
        const T solidity = net_solidity::Solidity(m_solidity, m_barD, m_barL0, frame.meshOpeningAngle);
        net_element_forces::AddNodeForces(m_panelForceParams, frame, law, fluid, solidity,
            velA_ned, velB_ned, velC_ned, waterVel_ned,
            nodeAforce_ned, nodeBforce_ned, nodeCforce_ned, addedLinearDrag, clampLog);
    }
 
    void NodeForceJacobian(const hydrodynamics::PanelLoadLaw& law, const hydrodynamics::Fluid& fluid,
        const double posA_ned[3], const double posB_ned[3], const double posC_ned[3],
        const double velA_ned[3], const double velB_ned[3], const double velC_ned[3],
        const double waterVel_ned[3], double addedLinearDrag, double dF[3][6][9]) const
    {
        using Dual = sfh::ad::Dual<12>;
        const int kVelocitySlot = 9;
 
        Dual pos[3][3], vel[3][3], water[3];
        const double* posIn[3] = {posA_ned, posB_ned, posC_ned};
        const double* velIn[3] = {velA_ned, velB_ned, velC_ned};
        for (int node = 0; node < 3; ++node) {
            for (int k = 0; k < 3; ++k) {
                pos[node][k] = Dual::Variable(posIn[node][k], 3 * node + k);
                vel[node][k] = Dual::Variable(velIn[node][k], kVelocitySlot + k);
            }
        }
        for (int k = 0; k < 3; ++k)
            water[k] = Dual(waterVel_ned[k]);
 
        Dual force[3][3];
        AddNodeForces<Dual>(law, fluid, pos[0], pos[1], pos[2], vel[0], vel[1], vel[2], water,
            force[0], force[1], force[2], addedLinearDrag);
 
        // CalcDampingForces applies the added drag only when it is positive.
        const double addedDragPerNode = addedLinearDrag > 0 ? addedLinearDrag / 3 : 0.0;
        const double damping          = m_panelForceParams.dampingCoeff;
        for (int node = 0; node < 3; ++node) {
            for (int j = 0; j < 3; ++j) {
                for (int k = 0; k < 3; ++k) {
                    for (int group = 0; group < 3; ++group)
                        dF[node][group][3 * j + k] = force[node][j].Grad(3 * group + k);
 
                    const double diagonal = (j == k) ? 1.0 : 0.0;
                    const double lawPart  = (force[node][j].Grad(kVelocitySlot + k) + diagonal * addedDragPerNode) / 3;
                    for (int other = 0; other < 3; ++other) {
                        const double self = (other == node) ? 1.0 : 0.0;
                        dF[node][3 + other][3 * j + k] = lawPart
                            + diagonal * (damping * (1.0 / 3.0 - self) - addedDragPerNode * self);
                    }
                }
            }
        }
    }
 
    hydrodynamics::PanelLoad<double> EvaluatePanelLoad(const hydrodynamics::PanelLoadLaw& law, const hydrodynamics::Fluid& fluid,
        const double* const posA_ned, const double* const posB_ned, const double* const posC_ned,
        const double elementVel[3], const double waterVel_ned[3], double& area) const
    {
        const net_element_forces::PanelFrame<double> frame =
            net_element_forces::ComputePanelFrame(m_panelForceParams, posA_ned, posB_ned, posC_ned);
        area = frame.area;
        return net_element_forces::EvaluatePanelLoad(m_panelForceParams, frame, law, fluid,
            Solidity(frame.meshOpeningAngle), elementVel, waterVel_ned);
    }
 
    bool HasBarDirections() const { return m_meshDet != 0; }
 
    double Solidity(double meshOpeningAngle) const { return net_solidity::Solidity(m_solidity, m_barD, m_barL0, meshOpeningAngle); }
 
    double Solidity(const double* const posA_ned, const double* const posB_ned, const double* const posC_ned) const
    {
        return Solidity(net_element_forces::ComputePanelFrame(m_panelForceParams, posA_ned, posB_ned, posC_ned).meshOpeningAngle);
    }
 
    double WakeSolidity(const double* const posA_ned, const double* const posB_ned, const double* const posC_ned) const
    {
        return net_solidity::WakeSolidity(m_solidity, m_barD, m_barL0,
                net_element_forces::ComputePanelFrame(m_panelForceParams, posA_ned, posB_ned, posC_ned).meshOpeningAngle);
    }
 
    double CalcNodeInertia(int node);
 
    double NumberOfTwineLengthsU() const;
 
    double TwineDiameter() const { return m_barD; }
 
#ifdef FH_VISUALIZATION
    void RenderInit(Ogre::Root* ogreRoot, int NumMeshPerLine = 1);
 
    void RenderUpdate( const double* const adPosA, const double* const adPosB, 
        const double* const adPosC );
#endif
 
 
private:
    void FillPanelForceParams();
 
    // Input 
    unsigned long  m_iD;                    
    std::string m_netName;              
    double m_nodeInertia[3];            
    //double m_rho_water;                   ///< Density of water.
    //double m_rho_twines;              ///< Material density of twines.
    //double m_EModulus;                    ///< E-modulus of twine material.
    double m_barD;                  
    double m_knotD;             
    double m_barL0;                 
    double m_nodeAPos_mesh[2];          
    double m_nodeBPos_mesh[2];          
    double m_nodeCPos_mesh[2];          
 
    // Calculated element properties
    //double m_twineCrossSecArea;           ///< The cross sectional area of each twine.
    double m_nodesDistAB_mesh[2];       
    double m_nodesDistAC_mesh[2];       
    double m_nodesDistBC_mesh[2];       
    double m_meshDet;                   
    double m_numMeshes;                 
    double m_numUBars;                  
    double m_numVBars;                  
    double m_numKnots;                  
    double m_weight;                    
    double m_eATwine;                   
 
    double m_NodeWeight;                
    double m_nodeAforce_ned[3];         
    double m_nodeBforce_ned[3];         
    double m_nodeCforce_ned[3];         
    double m_dampingCoeff;              
 
    NetPanelOptions m_options; 
    net_element_forces::PanelParams m_panelForceParams;
    net_solidity::PanelSolidity m_solidity;
#ifdef FILTER_FORCES
    double m_nodeAforce_ned_last[3];            
    double m_nodeBforce_ned_last[3];            
    double m_nodeCforce_ned_last[3];            
#endif
 
#ifdef FH_VISUALIZATION
    // Visualization specifics
    Ogre::Root*             m_ogreRoot;                 
    Ogre::SceneManager*     m_sceneMgr;                 
    Ogre::Entity*           m_renderElement[3];         
    Ogre::Entity*           m_renderSphere[3];
    Ogre::SceneNode*        m_renderNodeMesh;           
    Ogre::SceneNode*        m_renderNodeSphere[3];      
    Ogre::SceneNode*        m_renderNodeOutline;        
    Ogre::SceneNode*        m_renderNodeMaterial;       
#ifdef USE_DYNAMIC_LINES
    C3DLine*            m_lines = nullptr;                  
#else
    C3DLine*                m_lines = nullptr;                  
    Ogre::ColourValue       m_twineColour;              
#endif
    C3DLine*            m_lineOutline;              
    //Ogre::ManualObject *  m_netRenderObj;             ///< The manual object used for rendering the net element with a triangle skin.
 
    // Visualization variables
    unsigned short*         m_sortU = nullptr;                  
    unsigned short*         m_sortV = nullptr;                  
    unsigned long               m_numTwines;                
    unsigned long               m_numUTwines;               
    unsigned long               m_numVTwines;               
    double*                 m_interpCoeffAu = nullptr;          
    double*                 m_interpCoeffAv = nullptr;          
    double*                 m_interpCoeffBu = nullptr;          
    double*                 m_interpCoeffBv = nullptr;          
    int                     m_NumMeshPerLine;
#endif
 
};
 
 
#endif