DifferentialOperators.h

/*
 * Copyright 2025 Adobe. All rights reserved.
 * This file is licensed to you under the Apache License, Version 2.0 (the "License");
 * you may not use this file except in compliance with the License. You may obtain a copy
 * of the License at http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software distributed under
 * the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR REPRESENTATIONS
 * OF ANY KIND, either express or implied. See the License for the specific language
 * governing permissions and limitations under the License.
 */
#pragma once
 
#include <lagrange/SurfaceMesh.h>
 
#include <Eigen/Core>
#include <Eigen/Sparse>
 
namespace lagrange::polyddg {


template <typename Scalar, typename Index>
class DifferentialOperators
{
public:
    using Vector = Eigen::Matrix<Scalar, 1, 3>;
 
public:
    DifferentialOperators(SurfaceMesh<Scalar, Index>& mesh);
 
public:
    Eigen::SparseMatrix<Scalar> gradient() const;

    Eigen::SparseMatrix<Scalar> d0() const;

    Eigen::SparseMatrix<Scalar> d1() const;

    Eigen::SparseMatrix<Scalar> star0() const;

    Eigen::SparseMatrix<Scalar> star1(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> star2() const;

    Eigen::SparseMatrix<Scalar> flat() const;

    Eigen::SparseMatrix<Scalar> sharp() const;

    Eigen::SparseMatrix<Scalar> projection() const;

    Eigen::SparseMatrix<Scalar> inner_product_0_form() const;

    Eigen::SparseMatrix<Scalar> inner_product_1_form(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> inner_product_2_form() const;

    Eigen::SparseMatrix<Scalar> divergence(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> curl() const;

    Eigen::SparseMatrix<Scalar> laplacian(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> delta1(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> delta2() const;

    Eigen::SparseMatrix<Scalar> laplacian2() const;

    Eigen::SparseMatrix<Scalar> laplacian1(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> vertex_tangent_coordinates() const;

    Eigen::SparseMatrix<Scalar> facet_tangent_coordinates() const;

    Eigen::SparseMatrix<Scalar> levi_civita() const;

    Eigen::SparseMatrix<Scalar> levi_civita_nrosy(Index n) const;

    Eigen::SparseMatrix<Scalar> covariant_derivative() const;

    Eigen::SparseMatrix<Scalar> covariant_derivative_nrosy(Index n) const;

    Eigen::SparseMatrix<Scalar> shape_operator() const;

    Eigen::SparseMatrix<Scalar> adjoint_gradient() const;

    Eigen::SparseMatrix<Scalar> adjoint_shape_operator() const;

    Eigen::SparseMatrix<Scalar> connection_laplacian(Scalar lambda = 1) const;

    Eigen::SparseMatrix<Scalar> connection_laplacian_nrosy(Index n, Scalar lambda = 1) const;
 
public:
    Eigen::Matrix<Scalar, 3, Eigen::Dynamic> gradient(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> d0(Index fid) const;

    Eigen::Matrix<Scalar, 1, Eigen::Dynamic> d1(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, 3> flat(Index fid) const;

    Eigen::Matrix<Scalar, 3, Eigen::Dynamic> sharp(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> projection(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> inner_product_0_form(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> inner_product_1_form(
        Index fid,
        Scalar lambda = 1) const;

    Eigen::Matrix<Scalar, 1, 1> inner_product_2_form(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> laplacian(Index fid, Scalar lambda = 1)
        const;

    Eigen::Matrix<Scalar, 2, 2> levi_civita(Index fid, Index lv) const;

    Eigen::Matrix<Scalar, 2, 2> levi_civita_nrosy(Index fid, Index lv, Index n) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> levi_civita(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> levi_civita_nrosy(Index fid, Index n)
        const;

    Eigen::Matrix<Scalar, 4, Eigen::Dynamic> covariant_derivative(Index fid) const;

    Eigen::Matrix<Scalar, 4, Eigen::Dynamic> covariant_derivative_nrosy(Index fid, Index n) const;

    Eigen::Matrix<Scalar, 2, 2> shape_operator(Index fid) const;

    Eigen::Matrix<Scalar, 3, Eigen::Dynamic> adjoint_gradient(Index vid) const;

    Eigen::Matrix<Scalar, 2, 2> adjoint_shape_operator(Index vid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> covariant_projection(Index fid) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> covariant_projection_nrosy(
        Index fid,
        Index n) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> connection_laplacian(
        Index fid,
        Scalar lambda = 1) const;

    Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>
    connection_laplacian_nrosy(Index fid, Index n, Scalar lambda = 1) const;
 
private:
    void compute_vertex_normal_from_vector_area();
 
 
public:
    Eigen::Matrix<Scalar, 3, 2> facet_basis(Index fid) const;

    Eigen::Matrix<Scalar, 3, 2> vertex_basis(Index vid) const;
 
public:
    AttributeId get_vector_area_attribute_id() const { return m_vector_area_id; }

    AttributeId get_centroid_attribute_id() const { return m_centroid_id; }

    AttributeId get_vertex_normal_attribute_id() const { return m_vertex_normal_id; }
 
private:
    SurfaceMesh<Scalar, Index>& m_mesh;
    AttributeId m_vector_area_id = invalid_attribute_id();
    AttributeId m_centroid_id = invalid_attribute_id();
    AttributeId m_vertex_normal_id = invalid_attribute_id();
};

 
} // namespace lagrange::polyddg