Smoother.hpp

//---- GPL ---------------------------------------------------------------------
//
// Copyright (C)  Stichting Deltares, 2011-2021.
//
// 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 "MeshKernel/Constants.hpp"
 
namespace meshkernel
{
    class Mesh2D;
 
    class Smoother
    {
 
    public:
        explicit Smoother(const Mesh2D& mesh);
 
        void Compute();
 
        [[nodiscard]] auto GetWeight(UInt node, int connectedNode) const
        {
            return m_weights[node][connectedNode];
        }
 
        [[nodiscard]] auto GetConnectedNodeIndex(UInt node, int connectedNode) const
        {
            return m_connectedNodes[node][connectedNode];
        }
 
        [[nodiscard]] auto GetNumConnectedNodes(UInt node) const
        {
            return m_numConnectedNodes[node];
        }
 
    private:
        void Initialize();
 
        void ComputeTopologies();
 
        void ComputeOperators();
 
        bool ComputeCoordinates();
 
        void ComputeWeights();
 
        void ComputeOperatorsNode(UInt currentNode);
 
        void NodeAdministration(UInt currentNode);
 
        void ComputeNodeXiEta(UInt currentNode);
 
        [[nodiscard]] double OptimalEdgeAngle(UInt numFaceNodes,
                                              double theta1 = -1.0,
                                              double theta2 = -1.0,
                                              bool isBoundaryEdge = false) const;
 
        void AllocateNodeOperators(UInt topologyIndex);
 
        void SaveNodeTopologyIfNeeded(UInt currentNode);
 
        void ComputeJacobian(UInt currentNode, std::vector<double>& J) const;
 
        // The mesh to smooth
        const Mesh2D& m_mesh; 
 
        // Smoother weights
        std::vector<std::vector<double>> m_weights; 
 
        // Smoother operators
        std::vector<std::vector<std::vector<double>>> m_Gxi;  
        std::vector<std::vector<std::vector<double>>> m_Geta; 
        std::vector<std::vector<double>> m_Divxi;             
        std::vector<std::vector<double>> m_Diveta;            
        std::vector<std::vector<std::vector<double>>> m_Az;   
        std::vector<std::vector<double>> m_Jxi;               
        std::vector<std::vector<double>> m_Jeta;              
        std::vector<std::vector<double>> m_ww2;               
 
        // Smoother local caches
        std::vector<UInt> m_sharedFacesCache;                  
        std::vector<UInt> m_connectedNodesCache;               
        std::vector<std::vector<UInt>> m_faceNodeMappingCache; 
        std::vector<double> m_xiCache;                         
        std::vector<double> m_etaCache;                        
        std::vector<UInt> m_boundaryEdgesCache;                
        std::vector<double> m_leftXFaceCenterCache;            
        std::vector<double> m_leftYFaceCenterCache;            
        std::vector<double> m_rightXFaceCenterCache;           
        std::vector<double> m_rightYFaceCenterCache;           
        std::vector<double> m_xisCache;                        
        std::vector<double> m_etasCache;                       
 
        // Smoother topologies
        std::vector<UInt> m_nodeTopologyMapping;                               
        std::vector<std::vector<double>> m_topologyXi;                         
        std::vector<std::vector<double>> m_topologyEta;                        
        std::vector<std::vector<UInt>> m_topologySharedFaces;                  
        std::vector<std::vector<std::vector<UInt>>> m_topologyFaceNodeMapping; 
        std::vector<std::vector<UInt>> m_topologyConnectedNodes;               
 
        std::vector<UInt> m_numConnectedNodes;           
        std::vector<std::vector<UInt>> m_connectedNodes; 
 
        // Class variables
        UInt m_maximumNumConnectedNodes = 0; 
        UInt m_maximumNumSharedFaces = 0;    
 
        static constexpr int m_topologyInitialSize = 10; 
        static constexpr double m_thetaTolerance = 1e-4; 
    };
} // namespace meshkernel