SplineAlgorithms.hpp

//---- GPL ---------------------------------------------------------------------
//
// Copyright (C)  Stichting Deltares, 2011-2023.
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation version 3.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see <http://www.gnu.org/licenses/>.
//
// contact: delft3d.support@deltares.nl
// Stichting Deltares
// P.O. Box 177
// 2600 MH Delft, The Netherlands
//
// All indications and logos of, and references to, "Delft3D" and "Deltares"
// are registered trademarks of Stichting Deltares, and remain the property of
// Stichting Deltares. All rights reserved.
//
//------------------------------------------------------------------------------
 
#pragma once
 
#include <memory>
 
#include "MeshKernel/Constants.hpp"
#include "MeshKernel/Entities.hpp"
#include "MeshKernel/LandBoundary.hpp"
#include "MeshKernel/Utilities/LinearAlgebra.hpp"
 
namespace meshkernel
{
 
    class SplineAlgorithms
    {
    public:
        [[nodiscard]] static std::vector<Point> SecondOrderDerivative(const std::vector<Point>& splines, size_t startIndex, size_t endIndex);
 
        [[nodiscard]] static std::vector<double> SecondOrderDerivative(const std::vector<double>& coordinates, size_t startIndex, size_t endIndex);
 
        static std::tuple<Point, Point, double>
        ComputeCurvatureOnSplinePoint(const std::vector<Point>& splinePoints,
                                      const std::vector<Point>& splineDerivative,
                                      double adimensionalPointCoordinate,
                                      const Projection projection);
 
        static Point Evaluate(const std::vector<Point>& coordinates, const std::vector<Point>& secondDerivative, const double evaluationPoint);
 
        static void SnapSplineToBoundary(std::vector<Point>& splinePoints,
                                         const std::vector<Point>& splineDerivative,
                                         const LandBoundary& landBoundary,
                                         const Projection projection,
                                         const int numberOfIterations = constants::numeric::defaultSnappingIterations);
 
    private:
        static lin_alg::ColVector<double> ComputeSplineWeights(const lin_alg::ColVector<double>& xf,
                                                               const lin_alg::ColVector<double>& yf,
                                                               const Projection projection);
        static std::tuple<lin_alg::ColVector<double>, lin_alg::ColVector<double>> ComputeSamplePoints(const std::vector<Point>& splinePoints,
                                                                                                      const lin_alg::Matrix<double, Eigen::ColMajor>& aMatrix);
 
        static void SampleSpline(const std::vector<Point>& splinePoints,
                                 const size_t intermediatePointCount,
                                 std::vector<Point>& samplePoints);
 
        static void ComputeInterpolationMatrix(const Eigen::Index numberOfSplinePoints,
                                               const Eigen::Index intervalRefinement,
                                               Eigen::Index& numberOfSamplePoints,
                                               lin_alg::Matrix<double, Eigen::ColMajor>& interpolationMatrix);
 
        static lin_alg::Matrix<double, Eigen::ColMajor> ComputeLeastSquaresMatrixInverse(const lin_alg::Matrix<double, Eigen::ColMajor>& splineCoefficients,
                                                                                         const lin_alg::ColVector<double>& weights);
    };
 
} // namespace meshkernel