Cable.h

 
#ifndef CABLE_H
#define CABLE_H
#include <fhsim/simobject/StateFrame.h>
#include <string>
 
#include <fhsim/simobject/SimObject.h>
#include <fhsim_environment/ElementFlowInteraction.h>
#include <fhsim_environment/SeafloorForcesInterface.h>
#include <optional>
#include <vector>
namespace environment
{
class EnvironmentProvider;
}
class C3DLine;
 
class CCable : public SimObject
{
   public:
   CCable(const std::string& simObjectName, ISimObjectCreator* const creator);
 
   ~CCable(); 
 
   void OdeFcn(const double T, const double* const X, double* const XDot) const;
 
 
   void InitialConditionSetup(const double T, const double* const currentIC, double* const updatedIC, ISimObjectCreator* const creator);
 
   void FinalSetup(const double T, const double* const X, ISimObjectCreator* const creator); 
 
   void PreOdeFcn(const double T, const double* const X) override;
 
   const double* ForceA(const double T, const double* const X); 
   const double* ForceB(const double T, const double* const X); 
 
   void computeEndForces(const double T, const double* const X) const; 
 
   void computeEndpoints(const double T, const StateFrame& frame) const;
 
#ifdef FH_VISUALIZATION
   void RenderInit(Ogre::Root* const ogreRoot, ISimObjectCreator* const creator); 
   void RenderUpdate(const double T, const double* const X); 
#endif
 
   protected:
   void computeSegmentForces(const double T, int segment, const double* X, const double* pA, const double* pB, const double* vA, const double* vB, double* FA, double* FB, double waveBound, double* wetOut = nullptr, bool* inSeabedContactOut = nullptr) const;
   double nodeMass(double wetHalves) const;
   double waveInertiaCoefficient() const;
   std::optional<double> surfaceElevation(const double T, const double* pA, const double* pB, double waveBound) const;
   double wetFraction(const double* pA, const double* pB, std::optional<double> seaElevation, double* dWet_dpA = nullptr, double* dWet_dpB = nullptr) const;
 
 
   void computeVirtualPoint(const double* const target, const double* const point, const double* const targetVel, const double* const pointVel, double* virtualPoint, double* virtualPointVel, double residual, double dResidual) const;
 
   void testPositiveParam(double value, ISimObjectCreator* creator, std::string valueName); 
 
   static constexpr int kLiveSegment = -1;
 
   struct SegmentSample
   {
      double L = 0.0;
      double Q[3] = {0.0, 0.0, 0.0};
      double midPosition[3] = {0.0, 0.0, 0.0};
      double wet = 0.0; 
      double jacobianWet = 0.0; 
      double dWet_dpA[3] = {0.0, 0.0, 0.0}; 
      double dWet_dpB[3] = {0.0, 0.0, 0.0};
      bool hasAcceleration = false; 
      double acceleration[3] = {0.0, 0.0, 0.0};
      double depth = 0.0;
      double relativeHeight = 0.0;
      bool hasCurrent = false; 
      double seabedRamp = 0.0;
      double current[3] = {0.0, 0.0, 0.0};
   };
 
   SegmentSample sampleSegment(const double T, int segment, const double* pA, const double* pB, double waveBound, bool withWetDerivatives) const;
 
   void segmentForces(const SegmentSample& sample, const double* pA, const double* pB, const double* vA, const double* vB, double* FA, double* FB, double* wetOut = nullptr, bool* inSeabedContactOut = nullptr) const;
 
   struct SegmentJacobian
   {
      double dFA_dpA[9], dFA_dpB[9], dFA_dvA[9], dFA_dvB[9];
      double dFB_dpA[9], dFB_dpB[9], dFB_dvA[9], dFB_dvB[9];
   };
 
   void computeSegmentJacobian(const double T, int segment, const double* pA, const double* pB, const double* vA, const double* vB, SegmentJacobian& jac, double waveBound) const;
   void segmentJacobian(const SegmentSample& sample, const double* pA, const double* pB, const double* vA, const double* vB, SegmentJacobian& jac) const;
 
   void addSegmentHydrodynamicJacobian(const SegmentSample& sample, const double* vA, const double* vB, SegmentJacobian& jac) const;
 
   double m_cableLength; 
   double m_cableWeight; 
   double m_Emodulus; 
   double m_tensileCrossSectionArea; 
   double m_displacementCrossSectionArea; 
   double m_diameter; 
   double m_hydrodynamicDiameter; 
   double m_submergedWeight; 
   double m_waterDensity; 
   double m_internalDamping; 
   int m_numElements; 
   double m_elementMass; 
 
   double m_segmentMass; 
   double m_segmentSubmergedMass; 
   double m_segmentLength; 
   double m_segmentStiffness; 
   int* m_elementPositionIndex; 
   int* m_elementVelocityIndex; 
   int m_stateBase = 0; 
 
   ISignalPort* m_PositionA; 
   ISignalPort* m_PositionB; 
   ISignalPort* m_VelocityA; 
   ISignalPort* m_VelocityB; 
   ISignalPort* m_retractedLengthPortA; 
   ISignalPort* m_retractedLengthPortB; 
   ISignalPort* m_retractedSpeedPortA; 
   ISignalPort* m_retractedSpeedPortB; 
 
   mutable double m_endForceA[3]; 
   mutable double m_endForceB[3]; 
   mutable double m_nodeForceA[3]; 
   mutable double m_nodeForceB[3]; 
   mutable double m_wetEndA = 1.0; 
   mutable double m_wetEndB = 1.0; 
   mutable std::vector<double> m_nodeWetHalves; 
   mutable double m_retractedLengthA; 
   mutable double m_retractedLengthB; 
   mutable int m_retractedNodesA; 
   mutable int m_retractedNodesB; 
   mutable double m_virtualPointA[3]; 
   mutable double m_virtualPointB[3]; 
   mutable double m_virtualPointVelocityA[3]; 
   mutable double m_virtualPointVelocityB[3]; 
   mutable double m_residualA; 
   mutable double m_residualB; 
   ICommonComputation* m_computeEndForces; 
 
   mutable const double* m_posA; 
   mutable const double* m_posB; 
   mutable const double* m_velA; 
   mutable const double* m_velB; 
 
   environment::wake::ElementFlowInteraction m_flow;        
   environment::EnvironmentProvider*         m_environment; 
   environment::SeafloorForcesInterface*     m_seafloorForces = nullptr; 
   mutable bool m_endSegmentsInSeabedContact = false; 
   mutable bool m_inSeabedContact            = false; 
 
   public:
   bool HasJacobians() const override;
   void OdeJacobian(double T, const double* X, double* J, int nStates) override;
   int GetJacobianSparsity(int nStates, int* rowPtr, int* colIdx) override;
 
#ifdef FH_VISUALIZATION
   C3DLine* m_lines1 = nullptr; 
   C3DLine* m_lines2 = nullptr; 
#endif
};
#endif