From c6d1686c5d389e67871f72da428175be53895607 Mon Sep 17 00:00:00 2001 From: Emily Boudreaux Date: Mon, 9 Nov 2026 13:16:36 -0500 Subject: [PATCH] feat(mapping): The mapping code was all moved to the backend, this module now just brings the MFEM backend in The mapping code is fundamentally backend dependent. Previously we had DomainMapper here; however, since it depende on MFEM we have moved it to the backend. The mapping module now just includes parts of the backend so that prior code will not break --- src/include/serif/mapping/coefficients.hpp | 95 ------ .../compactification/compactification.hpp | 82 ----- .../mapping/compactification/input_result.hpp | 36 --- .../serif/mapping/compactification/kelvin.hpp | 197 ------------ .../mapping/compactification/variation.hpp | 35 --- src/include/serif/mapping/context.hpp | 92 ------ src/include/serif/mapping/data.hpp | 295 ------------------ src/include/serif/mapping/domain_mapper.hpp | 166 ---------- src/include/serif/mapping/evaluator.hpp | 93 ------ src/include/serif/mapping/mapping.hpp | 4 + src/include/serif/mapping/status.hpp | 58 ---- src/include/serif/mapping/transformations.hpp | 87 ------ src/include/serif/mapping/variation.hpp | 29 -- src/include/serif/mapping/volume_cache.hpp | 28 -- 14 files changed, 4 insertions(+), 1293 deletions(-) delete mode 100644 src/include/serif/mapping/coefficients.hpp delete mode 100644 src/include/serif/mapping/compactification/compactification.hpp delete mode 100644 src/include/serif/mapping/compactification/input_result.hpp delete mode 100644 src/include/serif/mapping/compactification/kelvin.hpp delete mode 100644 src/include/serif/mapping/compactification/variation.hpp delete mode 100644 src/include/serif/mapping/context.hpp delete mode 100644 src/include/serif/mapping/data.hpp delete mode 100644 src/include/serif/mapping/domain_mapper.hpp delete mode 100644 src/include/serif/mapping/evaluator.hpp create mode 100644 src/include/serif/mapping/mapping.hpp delete mode 100644 src/include/serif/mapping/status.hpp delete mode 100644 src/include/serif/mapping/transformations.hpp delete mode 100644 src/include/serif/mapping/variation.hpp delete mode 100644 src/include/serif/mapping/volume_cache.hpp diff --git a/src/include/serif/mapping/coefficients.hpp b/src/include/serif/mapping/coefficients.hpp deleted file mode 100644 index 204154a0..00000000 --- a/src/include/serif/mapping/coefficients.hpp +++ /dev/null @@ -1,95 +0,0 @@ -#pragma once - -/** - * @file coefficients.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the PhysicalPositionFunctionCoefficient class, which is a custom MFEM coefficient that evaluates a user-defined function at physical positions determined by a domain mapping. - * - * @todo This should probably be moved into a backend/mfem module - */ - -#include - -#include "mfem/fem/gridfunc.hpp" -#include "mfem/fem/coefficient.hpp" -#include "mfem/fem/eltrans.hpp" -#include "mfem/fem/intrules.hpp" -#include "mfem/linalg/vector.hpp" - -#include "serif/mapping/domain_mapper.hpp" -#include "serif/mapping/evaluator.hpp" - -namespace serif::mapping { - /** - * @class PhysicalPositionFunctionCoefficient - * - * @brief A custom MFEM coefficient that evaluates a user-defined function at physical positions determined by a domain mapping. - * - * @par Example Usage - * Lets say we have some analytic density profile of the form - * - * \f[ - * \rho(r) = \rho_0 \exp(-r^2 / R^2) - * \f] - * - * We might want to evaluate this function over some computational domain but in physical space since - * the analytic form is setup to be evaluated in physical space. We cannot simply calculate the - * position of each dof on the mesh and evaluate the function there since the mesh may be deformed. Rather - * if we use this coefficient we can offload the mapping to the mapping module - * - * \code{.cpp} - * // assuming there is already some defined mapper - * - * auto density = [](const mfem::Vector& x) -> double { - * double r = x.Norml2(); - * return rho0 * std::exp(-r*r / (R*R)); - * }; - * - * PhysicalPositionFunctionCoefficient physical_density(mapper, displacement, compactification_coordinate, density); - * // Then there would need to be some loopup code to get the rigth element transformation and integration point - * // for the physical position you want to evaluate the function at. - * double denity = physical_density.Eval(transformation, ip); - * - * \endcode - */ - class PhysicalPositionFunctionCoefficient : public mfem::Coefficient { - public: - using Func = std::function; ///< The type of user-defined function must be a scalar valued function of a 3D vector. - - /** - * @brief Construct a PhysicalPositionFunctionCoefficient with a given domain mapper, displacement, compactification coordinate, and user-defined function. - * - * @param mapper The domain mapper that defines the mapping from reference to physical coordinates. - * @param displacement The grid function representing the displacement field. - * @param compactification_coordinate The grid function representing the compactification coordinate. - * @param func The user-defined function to evaluate at physical positions. - * - * @todo The displacement grid function really should not be needed in this constructor since all required information - * is already in the mapper. Same goes for the compactification coordinate. I kept these to maintain compatibility - * with some older code; however, we should consider removing them. - */ - PhysicalPositionFunctionCoefficient( - const DomainMapper& mapper, - const mfem::GridFunction& displacement, - const mfem::GridFunction& compactification_coordinate, - Func func - ); - - /** - * @brief Evaluate the user-defined function at a given point in the reference domain, mapped to physical space. - * - * @param transformation The element transformation that defines the mapping from reference to physical coordinates. - * @param ip The integration point in the reference domain at which to evaluate the function. - * - * @return The value of the user-defined function evaluated at the physical position corresponding to the given reference point. - */ - double Eval( - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip - ) override; - private: - Func m_f; - GridFunctionMappingEvaluator m_mapping; ///< The mapping evaluator that computes physical positions from reference coordinates. - }; -} diff --git a/src/include/serif/mapping/compactification/compactification.hpp b/src/include/serif/mapping/compactification/compactification.hpp deleted file mode 100644 index d34851ec..00000000 --- a/src/include/serif/mapping/compactification/compactification.hpp +++ /dev/null @@ -1,82 +0,0 @@ -#pragma once - -/** - * @file compactification.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the definition of the ExteriorDomainMap interface, which is used to define the mapping from reference coordinates to physical coordinates in the exterior domain. - */ - -#include -#include - -#include "serif/mapping/status.hpp" -#include "serif/mapping/compactification/input_result.hpp" -#include "serif/mapping/compactification/variation.hpp" - - -namespace serif::mapping::compactification { - - /** - * @class ExteriorDomainMap - * - * @brief Abstract interface for some mapping from reference coordinates to physical coordinates in the exterior domain. - * Currently the only concrete implementation of this is the kelvin compactification. - */ - class ExteriorDomainMap { - public: - virtual ~ExteriorDomainMap() = default; - - /** - * @brief Check if the mapping is affine in the displacement, that is if the mapping is linear in the displacement and does not depend on the displacement. - * - * @return True if the mapping is affine in the displacement, false otherwise. - */ - [[nodiscard]] virtual bool is_affine_in_displacement() const noexcept = 0; - - /** - * @brief Evaluate the mapping at a given point in the reference domain. - * - * @param input The input to the mapping, which contains the reference coordinates and the displacement. - * @param result The result of the mapping, which contains the physical coordinates and the Jacobian of the mapping. - * - * @return The status of the mapping evaluation, which indicates if the evaluation was successful or if there was an error. - * - */ - [[nodiscard]] virtual MappingStatus evaluate( - const ExteriorMapInput &input, - ExteriorMapResult &result - ) const = 0; - - /** - * @brief Evaluate the variation of the mapping at a given point in the reference domain. - * - * @param input The input to the mapping, which contains the reference coordinates and the displacement. - * @param result The result of the mapping, which contains the physical coordinates and the Jacobian of the mapping. - * @param direction The direction of the variation, which is a vector in the reference domain. - * @param variation The variation of the mapping, which contains the variation of the physical coordinates and the variation of the Jacobian. - * - * @return The status of the mapping evaluation, which indicates if the evaluation was successful or if there was an error. - */ - [[nodiscard]] virtual MappingStatus evaluate_variation( - const ExteriorMapInput &input, - const ExteriorMapResult &result, - const ExteriorMapDirection &direction, - ExteriorMapVariation &variation - ) const = 0; - - /** - * @brief Return the name of the mapping, which is used for debugging and logging. - * - * @return The name of the mapping. - */ - [[nodiscard]] virtual std::string_view get_name() const noexcept = 0; - - /** - * @brief Clone the mapping, which creates a new instance of the mapping with the same parameters. - * - * @return A unique pointer to the cloned mapping. - */ - [[nodiscard]] virtual std::unique_ptr clone() const = 0; - }; -} \ No newline at end of file diff --git a/src/include/serif/mapping/compactification/input_result.hpp b/src/include/serif/mapping/compactification/input_result.hpp deleted file mode 100644 index 561c71a7..00000000 --- a/src/include/serif/mapping/compactification/input_result.hpp +++ /dev/null @@ -1,36 +0,0 @@ -#pragma once - -/** - * @file input_result.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the input and result structures for the compactification mapping. - */ - -#include "mfem/linalg/vector.hpp" -#include "mfem/linalg/densemat.hpp" - -namespace serif::mapping::compactification { - /** - * @struct ExteriorMapInput - * - * @brief Structure representing the input to the exterior domain mapping. This includes the reference position, displaced position, displacement Jacobian, compactification coordinate gradient, and compactification coordinate. - */ - struct ExteriorMapInput { - const mfem::Vector &reference_position; ///< Position in reference space to evaluate the mapping at. - const mfem::Vector &displaced_position; ///< Position in displaced space to evaluate the mapping at (note this is displaced but not yet compactified). - const mfem::DenseMatrix &displacement_jacobian; ///< Jacobian of the displacement mapping at the reference position. - const mfem::Vector &compactification_coordinate_gradient; ///< Gradient of the compactification coordinate at the reference position. - double compactification_coordinate; ///< The compactification coordinate at the reference position. This ranges from 0 to 1, where 0 is the inner surface and 1 is the outer surface of the reference domain. - }; - - /** - * @struct ExteriorMapResult - * - * @brief Structure representing the result of the exterior domain mapping. This includes the physical position and the Jacobian of the mapping. - */ - struct ExteriorMapResult { - mfem::Vector physical_position; ///< Position in physical space after applying the mapping. - mfem::DenseMatrix mapping_jacobian; ///< Jacobian of the mapping at the reference position, which maps from reference coordinates to physical coordinates. - }; -} diff --git a/src/include/serif/mapping/compactification/kelvin.hpp b/src/include/serif/mapping/compactification/kelvin.hpp deleted file mode 100644 index 40dd68e8..00000000 --- a/src/include/serif/mapping/compactification/kelvin.hpp +++ /dev/null @@ -1,197 +0,0 @@ -#pragma once - -/** - * @file kelvin.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the definition of the KelvinCompactification class, which is a concrete implementation of the ExteriorDomainMap interface for the Kelvin compactification. - */ - -#include - -#include "mfem/linalg/vector.hpp" - -#include "serif/mapping/status.hpp" -#include "serif/mapping/compactification/compactification.hpp" -#include "serif/utils/user_options/mapping/compactification/kelvin.hpp" - -namespace serif::mapping::compactification { - - /** - * @class KelvinCompactification - * - * @brief Concrete implementation of the ExteriorDomainMap interface for the Kelvin compactification. - * - * Broadly, for some reference position \f$\mathbf{x_{r}}\f$ on some reference domain \f$\Omega_{r}\f$ the Kelvin compactification maps to a physical position \f$\mathbf{x_{p}}\f$ on some physical domain \f$\Omega_{p}\f$ through the mapping - * \f[ - * \mathbf{x_{p}} = \frac{\mathbf{x_{r}}}{1 - \frac{|\mathbf{x_{r}}|}{R_{\infty}}} - * \f] - * - * This results in the inner surface of the domain being mapped to the exterior surface of the physical domain while the outer - * surface of the exterior domain is mapped to infinity. There are some implications of this, namely that we - * loose accuracy in all things as we approach the outer surface of the domain. However, it is nice in so far - * as we can apply true boundary conditions at infinity. - */ - class KelvinCompactification final : public ExteriorDomainMap { - public: - /** - * @brief Construct a KelvinCompactification with the given options. - * - * @param options The options for the Kelvin compactification, which must be valid. - * - * @throws std::invalid_argument If the reference stellar radius is not finite or positive - * @throws std::invalid_argument If the reference infinity radius is not finite or positive - * @throws std::invalid_argument If the reference infinity radius is not greater than the reference stellar radius - * @throws std::invalid_argument If the coordinate tolerance is not finite or positive or is greater than or equal to one - */ - explicit KelvinCompactification( - const utils::user_options::KelvinCompactificationOptions& options - ); - - /** - * @brief Check if the mapping is affine in the displacement. For the Kelvin compactification, this is always true - * - * @return false - */ - [[nodiscard]] bool is_affine_in_displacement() const noexcept override; - - /** - * @brief Evaluate the mapping at a given point in the reference domain. - * - * @param[in] input The input point in the reference domain to evaluate the mapping at. - * @param[out] result The result of the mapping evaluation, which will be filled in by the method. - * - * @return A status indicating the success or failure of the mapping evaluation. Mapping may fail if - * - the size of the input reference position is less than or equal to zero - * - the input displaced position size does not match the input reference position size - * - the input compactification coordinate gradient size does not match the input reference position size - * - the displacement_jacobian height does not match the input reference position size - * - the displacement_jacobian width does not match the input reference position size - * - the input reference, displaced, compactification coordinate gradient or jacobian contains non-finite values - * - There is some issue computing the radial factors (this error will be passed through) - * - * The mapping jacobian comes from trivial application of the product rule - * recall: \f$r_i = scale * x_i\f$ where \f$r\f$ is the physical position and \f$x\f$ is the displaced position. - * then we can differentiate wrt. \f$X_j\f$ holding nothing fixed. Note the capital \f$X\f$ here, this is the mesh coordinate not the displaced position. - * Lets call this jacobian \f$F\f$ - * \f[ - * F = \frac{\partial r_i}{\partial X_{j}} - * \f] - * - * \f$F\f$ then tells us how the physical position changes as we move along mesh coordinates - * - * Lets then apply this to the function we have for the kelvin compactification - * \f[ - * F = scale * \frac{\partial x_i}{\partial X_j} + x_i * \frac{\partial scale}{\partial X_j} - * \f] - * - * In the implementation you can see the displacement jacobian (\f$\frac{\partial x_i}{\partial X_j}\f$) and the scale derivative (\f$\frac{\partial scale}{\partial X_j}\f$) being applied to compute the mapping jacobian. - */ - [[nodiscard]] MappingStatus evaluate( - const ExteriorMapInput& input, - ExteriorMapResult& result - ) const override; - - - /** - * @brief Evaluate the variation of the mapping at a given point in the reference domain. - * - * @param[in] input The input point in the reference domain to evaluate the mapping at. - * @param[in] result The result of the mapping evaluation, which will be used in the variation computation. - * @param[in] direction The direction in which to evaluate the variation. - * @param[out] variation The result of the variation evaluation, which will be filled in by the method. - * - * @return A status indicating the success or failure of the variation evaluation. The reasons for a failure are effectively the same - * as the reasons for a failure in the evaluate method just extended to cover the additional inputs. - */ - [[nodiscard]] MappingStatus evaluate_variation( - const ExteriorMapInput& input, - const ExteriorMapResult& result, - const ExteriorMapDirection& direction, - ExteriorMapVariation& variation - ) const override; - - /** - * @brief Get the name of the mapping, which is "KelvinCompactification" for this class. - * @return The name of the mapping. - */ - [[nodiscard]] std::string_view get_name() const noexcept override; - - /** - * @brief Clone the mapping, returning a unique pointer to a new instance of the mapping. - * @return A unique pointer to a new instance of the mapping. - */ - [[nodiscard]] std::unique_ptr clone() const override; - - /** - * @brief Get the reference stellar radius for the mapping, which is the radius of the inner surface of the reference domain. - * @return The reference stellar radius for the mapping. - */ - [[nodiscard]] double get_reference_stellar_radius() const noexcept; - - /** - * @brief Get the reference infinity radius for the mapping, which is the radius of the outer surface of the reference domain. - * @return The reference infinity radius for the mapping. - */ - [[nodiscard]] double get_reference_infinity_radius() const noexcept; - - /** - * @brief Get the coordinate tolerance for the mapping, which is the tolerance for the compactification coordinate. - * @return The coordinate tolerance for the mapping. - */ - [[nodiscard]] double get_coordinate_tolerance() const noexcept; - - /** - * @brief Evaluate the mapping jacobian variation factors at some given compactification coordinate. - * - * @param[in] compactification_coordinate The compactification coordinate at which to evaluate the mapping jacobian variation factors. (Note this runs from 0 to 1, where 0 is the inner surface and 1 is the outer surface of the reference domain) - * @param[in] coordinateGradient The gradient of the compactification coordinate with respect to the mesh coordinates, which must be finite and have size 3. - * - * @param[out] scale The scale factor for the mapping jacobian, which will be filled in by the method. - * @param[out] scaleGradient The gradient of the scale factor with respect to the mesh coordinates, which will be filled in by the method. - * - * @return The status of the mapping evaluation. - */ - [[nodiscard]] MappingStatus evaluate_jacobian_variation_factors( - double compactification_coordinate, - const mfem::Vector &coordinateGradient, - double &scale, - mfem::Vector &scaleGradient - ) const; - - private: - /** - * @struct RadialFactors - * - * @brief A struct to hold the radial factors for the Kelvin compactification mapping. - */ - struct RadialFactors { - double compactification_coordinate; ///< The compactification coordinate, which runs from 0 to 1, where 0 is the inner surface and 1 is the outer surface of the reference domain. - double computational_radius; ///< The computational radius, which is the radius of the reference domain at the given compactification coordinate. - double scale; ///< The scale factor for the mapping jacobian. - double scale_derivative; ///< The derivative of the scale factor with respect to the compactification coordinate. - }; - - private: - /** - * @brief Compute the radial factors for the Kelvin compactification mapping at a given compactification coordinate. - * - * @param[in] compactification_coordinate The compactification coordinate at which to compute the radial factors. (Note this runs from 0 to 1, where 0 is the inner surface and 1 is the outer surface of the reference domain) - * @param[out] radial_factors The radial factors, which will be filled in by the method. - * - * @return The status of the mapping evaluation. - * - * @note See the implementation file for more details on the actual algorithm used here. - */ - [[nodiscard]] MappingStatus compute_radial_factors( - double compactification_coordinate, - RadialFactors &radial_factors - ) const; - - private: - utils::user_options::KelvinCompactificationOptions m_options; ///< The options for the Kelvin compactification, which must be valid. - const std::string_view m_name = "KelvinCompactification"; ///< The name of the mapping, which is "KelvinCompactification" for this class. - }; - - -} \ No newline at end of file diff --git a/src/include/serif/mapping/compactification/variation.hpp b/src/include/serif/mapping/compactification/variation.hpp deleted file mode 100644 index 96196c28..00000000 --- a/src/include/serif/mapping/compactification/variation.hpp +++ /dev/null @@ -1,35 +0,0 @@ -#pragma once - -/** - * @file variation.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the structures for the variation of the compactification mapping. - */ - -#include "mfem/linalg/vector.hpp" -#include "mfem/linalg/densemat.hpp" - -namespace serif::mapping::compactification { - /** - * @struct ExteriorMapDirection - * - * @brief Structure representing the direction of the variation in the exterior domain mapping. - * This includes the variation of the displaced position and the variation of the displacement Jacobian. - */ - struct ExteriorMapDirection { - const mfem::Vector &displaced_position_variation; ///< Variation of the position in displaced space. - const mfem::DenseMatrix &displacement_jacobian_variation; ///< Variation of the Jacobian of the displacement mapping. - }; - - /** - * @struct ExteriorMapVariation - * - * @brief Structure representing the variation of the exterior domain mapping. - * This includes the variation of the physical position and the variation of the mapping Jacobian. - */ - struct ExteriorMapVariation { - mfem::Vector physical_position_variation; ///< Variation of the position in physical space. - mfem::DenseMatrix mapping_jacobian_variation; ///< Variation of the Jacobian of the mapping. - }; -} diff --git a/src/include/serif/mapping/context.hpp b/src/include/serif/mapping/context.hpp deleted file mode 100644 index 6d8f3f3e..00000000 --- a/src/include/serif/mapping/context.hpp +++ /dev/null @@ -1,92 +0,0 @@ -#pragma once - -/** - * @file context.hpp - * @author Emily M. Boudreaux - * - * @brief File containing the mapping context structures for volume and face mappings. - * - * @todo Like a lot of the parts of the mapping module written before I started factoring mfem stuff into - * backend modules, a good amount of this might be movable to some mfem backend - * - * All of these different "context" structures are built to contain sufficient information that the - * caller which gets one of these back can evaluate the mapping and its derivatives properly - * at some point. - */ - -#include "mfem/linalg/densemat.hpp" -#include "mfem/linalg/vector.hpp" - - -namespace serif::mapping { - /** - * @struct VolumeQuadratureContext - * - * @brief A structure representing the context of a volume quadrature point, including the inverse Jacobian of the mapping and the determinant of the mapping. - * - * @todo Add a worked example here - */ - struct VolumeQuadratureContext { - mfem::DenseMatrix J_element_to_physical_inverse; ///< The inverse Jacobian of the mapping from the element to the physical space. - double mapping_determinate; ///< The determinant of the mapping from the element to the physical space, - double physical_integration_weight; ///< The physical integration weight for the quadrature point. - }; - - /** - * @struct FaceQuadratureContext - * - * @brief A structure representing the context of a face quadrature point, including the normal vector, surface weight, and scale factor for the dot product of the velocity and normal vector. - * - * @todo Add a worked example here - */ - struct FaceQuadratureContext { - mfem::Vector normal; ///< The normal vector at the face quadrature point. - double ds; ///< The surface weight for the face quadrature point. - double v_dot_n_scale; ///< The scale factor for the dot product of the velocity and normal vector at the face quadrature point. - }; - - /** - * @struct MappingPointContext - * - * @brief A structure representing the context of a mapping point, including the reference position, displaced position, physical position, displacement Jacobian, mapping Jacobian, inverse mapping Jacobian, mapping determinant, and compactification status. - * - * @todo Add a worked example here - */ - struct MappingPointContext { - mfem::Vector reference_position; ///< The reference position of the mapping point. - mfem::Vector displaced_position; ///< The displaced position of the mapping point. - mfem::Vector physical_position; ///< The physical position of the mapping point. - mfem::DenseMatrix displacement_jacobian; ///< The displacement Jacobian of the mapping point. - mfem::DenseMatrix mapping_jacobian; ///< The mapping Jacobian of the mapping point. - mfem::DenseMatrix inverse_mapping_jacobian; ///< The inverse mapping Jacobian of the mapping point. - double mapping_determinate{0.0}; ///< The determinant of the mapping Jacobian at the mapping point. - bool compactified{false}; ///< A flag indicating whether the mapping point is compactified (true) or not (false). - }; - - /** - * @struct VolumeMappingContext - * - * @brief A structure representing the context of a volume mapping, including the mapping point context and volume quadrature context. - * - * @todo Add a worked example here - */ - struct VolumeMappingContext { - MappingPointContext mapping; ///< The mapping point context for the volume mapping. - VolumeQuadratureContext quadrature; ///< The volume quadrature context for the volume mapping. - }; - - /** - * @struct FaceMappingContext - * - * @brief A structure representing the context of a face mapping, including the mapping point context and face quadrature context. - * - * @todo Add a worked example here - */ - struct FaceMappingContext { - MappingPointContext mapping; ///< The mapping point context for the face mapping. - FaceQuadratureContext quadrature; ///< The face quadrature context for the face mapping. - mfem::Vector reference_normal; ///< The reference normal vector for the face mapping. - double reference_surface_weight{0.0}; ///< The reference surface weight for the face mapping. - double physical_surface_weight{0.0}; ///< The physical surface weight for the face mapping. - }; -} \ No newline at end of file diff --git a/src/include/serif/mapping/data.hpp b/src/include/serif/mapping/data.hpp deleted file mode 100644 index 8c378f37..00000000 --- a/src/include/serif/mapping/data.hpp +++ /dev/null @@ -1,295 +0,0 @@ -#pragma once - -/** - * @file data.hpp - * @author Emily M. Boudreaux - * - * @brief File containing some misc data structures for the mapping module. This was the very first file - * I ported over from MeanField and I think because of that its organization is actually quite wonky. - * For one its int this poorly named "data" file, which could mean any number of things. Regardless, - * these are basically low leval containers for useful stuff and also some scratch workspaces - * for mapping operations. - * - * @todo As said above there is some refactoring to make this cleaner which could happen. - */ - -#include "mfem/fem/fe/fe_base.hpp" -#include "mfem/fem/fespace.hpp" -#include "mfem/linalg/vector.hpp" -#include "mfem/linalg/densemat.hpp" -#include "mfem/mesh/mesh.hpp" - -#include -#include -#include - -#include "serif/mapping/compactification/kelvin.hpp" -#include "serif/mapping/context.hpp" -#include "serif/mapping/data.hpp" - - -namespace serif::mapping { - /** - * @enum FaceElementSide - * - * @brief An enumeration representing the side of a face element, either ELEMENT_1 or ELEMENT_2. There is no inherent ordering - * to these; however, broadly I have been thinking about ELEMENT_1 as the "inner" side and ELEMENT_2 as the "outer" side. - */ - enum class FaceElementSide : uint8_t { - ELEMENT_1, ///< Represents the first side of a face element. - ELEMENT_2 ///< Represents the second side of a face element. - }; - - /** - * @class ElementDisplacementData - * - * @brief A class representing the displacement data for a finite element, including the finite element, degree of freedom matrix, dimension, and ordering. - */ - class ElementDisplacementData { - public: - - /** - * @brief Constructs an ElementDisplacementData object with the given finite element, displacement degrees of freedom, and ordering. - * - * @param element The finite element associated with the displacement data. - * @param displacement_dofs The vector containing the displacement degrees of freedom. - * @param ordering The ordering type for the degrees of freedom (default is byNODES). - * - * @throws std::invalid_argument if the finite element has zero or negative degrees of freedom, - * or if the size of the displacement_dofs vector is not a positive multiple of the element's degree of freedom count. - * - * If ordering is set to byNODES then the dof matrix will be populated such that for row i and column j, with n total dofs, the value is set to - * dof_matrix(i, j) = displacement_dofs(i + j * n). Otherwise dof_matrix(i, j) = displacement_dofs(j + i * dim) for some number of dimensions dim (read as the size of the displacement_dof vector / dof count - * and implicitly truncated to an integer). - */ - ElementDisplacementData( - const mfem::FiniteElement& element, - mfem::Vector& displacement_dofs, - mfem::Ordering::Type ordering = mfem::Ordering::byNODES - ); - - /** - * @brief Returns the finite element associated with the displacement data. - * @return A constant reference to the finite element. - */ - [[nodiscard]] const mfem::FiniteElement& get_element() const noexcept; - - /** - * @brief Returns the degree of freedom matrix for the displacement data. - * @return A constant reference to the degree of freedom matrix. - */ - [[nodiscard]] const mfem::DenseMatrix &get_dof_matrix() const noexcept; - - /** - * @brief Returns the dimension of the displacement data. - * @return The dimension of the displacement data. - */ - [[nodiscard]] int get_dimension() const noexcept; - - /** - * @brief Returns the degree of freedom count for the displacement data. - * @return The degree of freedom count for the displacement data. - */ - [[nodiscard]] int get_dof_count() const noexcept; - - /** - * @brief Returns the ordering type for the degrees of freedom in the displacement data. - * @return The ordering type for the degrees of freedom. - */ - [[nodiscard]] mfem::Ordering::Type get_ordering() const noexcept; - - private: - const mfem::FiniteElement *m_element; ///< Pointer to the finite element associated with the displacement data. - mfem::DenseMatrix m_dof_matrix; ///< Matrix representing the degrees of freedom for the displacement data. - int m_dimension; ///< The dimension of the displacement data (number of components). - mfem::Ordering::Type m_ordering; ///< The ordering type for the degrees of freedom (byNODES or byVDIM). - }; - - /** - * @struct CompactificationPointData - * - * @brief A structure representing the compactification point data, including the compactification coordinate and coordinate gradient. - */ - struct CompactificationPointData { - double compactification_coordinate{0.0}; ///< The compactification coordinate value at the point. - mfem::Vector coordinate_gradient; ///< The gradient of the compactification coordinate at the point. - }; - - /** - * @brief Creates an ElementDisplacementData object from a finite element and a vector of displacement degrees of freedom. - * - * @param element The finite element associated with the displacement data. - * @param displacement_dofs The vector containing the displacement degrees of freedom. - * @return An ElementDisplacementData object constructed from the provided finite element and displacement degrees of freedom. - * - * @note This is just a convenience function which calls the constructor for ElementDisplacementData with the ordering set to byNODES. If you want to use a different ordering, call the constructor directly. - */ - [[nodiscard]] ElementDisplacementData element_displacement_data_from_element_vdofs( - const mfem::FiniteElement& element, - mfem::Vector &displacement_dofs - ); - - /** - * @class ElementCompactificationData - * - * @brief A class representing the compactification data for a finite element, including the finite element and degrees of freedom. - */ - class ElementCompactificationData { - public: - /** - * @brief Constructs an ElementCompactificationData object with the given finite element and degrees of freedom. - * - * @param element The finite element associated with the compactification data. - * @param dofs The vector containing the degrees of freedom for the compactification data. - * - * @throws std::invalid_argument if the finite element is not a scalar finite element, does not provide a gradient, or has zero or negative degrees of freedom, - * or if the size of the dofs vector does not match the finite element's degree of freedom count. - */ - ElementCompactificationData( - const mfem::FiniteElement& element, - const mfem::Vector& dofs - ); - - /** - * @brief Returns the finite element associated with the compactification data. - * @return A constant reference to the finite element. - */ - [[nodiscard]] const mfem::FiniteElement& get_element() const noexcept; - - /** - * @brief Returns the degrees of freedom for the compactification data. - * @return A constant reference to the vector of degrees of freedom. - */ - [[nodiscard]] const mfem::Vector &get_dofs() const noexcept; - - /** - * @brief Returns the degree of freedom count for the compactification data. - * @return The degree of freedom count for the compactification data. - */ - [[nodiscard]] int get_dof_count() const noexcept; - - private: - const mfem::FiniteElement *m_element; ///< Pointer to the finite element associated with the compactification data. - mfem::Vector m_dofs; ///< Vector representing the degrees of freedom for the compactification data. - }; - - struct ElementMappingData { - const ElementDisplacementData &displacement; - const ElementCompactificationData &compactification; - }; - - /** - * @class DMWorkspace - * - * @brief A class representing a workspace for displacement mapping operations, including various vectors and matrices used in the mapping process. - * primarily this is just a separate class so we can more easily make this code thread-safe in the future. At construction of - * an instance of this class loads of Heap allocation is done (or really whatever mfem does internally for allocation for vector which dependents on - * the selected device) and then that allocation can be reused. - */ - class DMWorkspace { - public: - /** - * @brief Constructs a DMWorkspace object with the specified dimension. - * - * @param dimension The dimension of the workspace (default is 3). - * - * @throws std::invalid_argument if the specified dimension is not positive and finite. - */ - explicit DMWorkspace (int dimension = 3); - - /** - * @brief Sets the dimension of the workspace. - * - * @param dimension The new dimension of the workspace. - * - * @throws std::invalid_argument if the specified dimension is not positive and finite. - * - * @post All internal vectors are sized to dimension - * @post All internal matrices are sized to dimension x dimension - * @post The compactification coordinate is set to 0.0 - */ - void set_dimension(int dimension); - - /** - * @brief Returns the dimension of the workspace. - * @return The dimension of the workspace. - */ - [[nodiscard]] int get_dimension() const noexcept; - - public: // This is a workspace class, these members are kept public so we do not need indirection through getters to access them - mfem::Vector m_shape; ///< Pre-allocated vector for shape function values at a quadrature point. - mfem::DenseMatrix m_reference_dshape; ///< Pre-allocated matrix for reference shape function derivatives at a quadrature point. - mfem::DenseMatrix m_mesh_dshape; ///< Pre-allocated matrix for mesh shape function derivatives at a quadrature point. - mfem::DenseMatrix m_reference_field_jacobian; ///< Pre-allocated matrix for the reference field Jacobian at a quadrature point. - mfem::Vector m_field_value; ///< Pre-allocated vector for field values at a quadrature point. - mfem::DenseMatrix m_field_jacobian; ///< Pre-allocated matrix for field Jacobian at a quadrature point. - - mfem::Vector m_compactification_shape; ///< Pre-allocated vector for compactification shape function values at a quadrature point. - mfem::DenseMatrix m_compactification_dshape; ///< Pre-allocated matrix for compactification shape function derivatives at a quadrature point. - CompactificationPointData m_compactification_point; ///< Pre-allocated data for compactification point. - - mfem::Vector m_reference_normal; ///< Pre-allocated vector for reference normal at a quadrature point. - mfem::Vector m_mapped_normal; ///< Pre-allocated vector for mapped normal at a quadrature point. - mfem::DenseMatrix m_full_element_jacobian; ///< Pre-allocated matrix for full element Jacobian at a quadrature point. - - mfem::Vector m_vector_temp; ///< Pre-allocated temporary vector for intermediate calculations. - mfem::DenseMatrix m_matrix_temp_1; ///< Pre-allocated temporary matrix for intermediate calculations. - mfem::DenseMatrix m_matrix_temp_2; ///< Pre-allocated temporary matrix for intermediate calculations. - - compactification::ExteriorMapResult m_exterior_result; ///< Pre-allocated result for exterior mapping at a quadrature point. - compactification::ExteriorMapVariation m_exterior_variation; ///< Pre-allocated variation for exterior mapping at a quadrature point. - private: - int m_dimension; ///< The number of spatial dimensions for the workspace (e.g., 2D or 3D). - }; - - /** - * @struct DMFixedMeshGeometry - * - * @brief A structure representing the fixed mesh geometry at a specific point, including position, Jacobian, inverse Jacobian, and weight. - */ - struct DMFixedMeshGeometry { - mfem::Vector position; ///< The physical position of the point in the mesh. - mfem::DenseMatrix jacobian; ///< The Jacobian matrix of the mapping from the reference element to the physical element at the point. - mfem::DenseMatrix inverse_jacobian; ///< The inverse of the Jacobian matrix of the mapping from the reference element to the physical element at the point. - double weight{0.0}; ///< The mapping weight associated with the point, typically used for integration purposes. - }; - - /** - * @struct DMGeometryCacheStatistics - * - * @brief A structure representing statistics for the geometry cache, including hit and miss counts and the number of points cached. - */ - struct DMGeometryCacheStatistics { - std::uint64_t hits{0}; ///< The number of cache hits, indicating how many times a requested geometry was found in the cache. - std::uint64_t misses{0}; ///< The number of cache misses, indicating how many times a requested geometry was not found in the cache and had to be computed. - std::size_t points{0}; ///< The number of points currently stored in the cache. - }; - - /** - * @struct DMReferenceDisplacementSample - * - * @brief A structure representing a sample of the reference displacement, including the value and gradient of the displacement at a specific point. - */ - struct DMReferenceDisplacementSample { - const mfem::Vector& value; ///< The value of the reference displacement at the point. - const mfem::DenseMatrix& gradient; ///< The gradient of the reference displacement at the point. - }; - - /** - * @struct DMGeometryCache - * - * @brief A structure representing a cache for geometry data, including the mesh, sequence numbers, and cached points. - */ - struct DMGeometryCache { - bool enabled{true}; ///< A flag indicating whether the geometry cache is enabled (true) or disabled (false). - const mfem::Mesh* mesh{nullptr}; ///< A pointer to the mesh associated with the cached geometry data. - - long sequence{-1}; ///< The sequence number of the cached geometry data. - long nodesSequence{-1}; ///< The sequence number of the cached nodes data. - - std::map, DMFixedMeshGeometry> points; ///< A map storing the cached geometry points. - - uint64_t hits{0}; ///< The number of cache hits. - uint64_t misses{0}; ///< The number of cache misses. - }; -} diff --git a/src/include/serif/mapping/domain_mapper.hpp b/src/include/serif/mapping/domain_mapper.hpp deleted file mode 100644 index ad9ad703..00000000 --- a/src/include/serif/mapping/domain_mapper.hpp +++ /dev/null @@ -1,166 +0,0 @@ -#pragma once - -#include -#include - -#include "mfem/linalg/vector.hpp" -#include "mfem/linalg/densemat.hpp" -#include "mfem/fem/eltrans.hpp" -#include "mfem/fem/intrules.hpp" - -#include "serif/mapping/compactification/kelvin.hpp" -#include "serif/mapping/context.hpp" -#include "serif/mapping/data.hpp" -#include "serif/mapping/variation.hpp" - -#include "serif/utils/user_options/mapping/domain_mapper.hpp" - - -namespace serif::mapping { - class DomainMapper { - public: - DomainMapper(const DomainMapper&) = delete; - DomainMapper& operator=(const DomainMapper&) = delete; - - ~DomainMapper() = default; - DomainMapper(DomainMapper&&) = default; - DomainMapper& operator=(DomainMapper&&) = default; - - - DomainMapper( - utils::user_options::DomainMapperOptions options, - std::unique_ptr exterior_map - ); - - std::unique_ptr clone() const; - - void set_geometry_cache_enabled(bool enabled) const; - [[nodiscard]] bool is_geometry_cache_enabled() const noexcept; - void clear_geometry_cache() const; - [[nodiscard]] DMGeometryCacheStatistics get_geometry_cache_statistics() const; - - [[nodiscard]] const DMFixedMeshGeometry* get_fixed_mesh_geometry(mfem::ElementTransformation& transformation, const mfem::IntegrationPoint& ip) const; - - [[nodiscard]] MappingStatus evaluate_point( - const ElementMappingData& element_data, - mfem::ElementTransformation &transformation, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - MappingPointContext& context - ) const; - - [[nodiscard]] MappingStatus evaluate_volume( - const ElementMappingData& element_data, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - VolumeMappingContext& context, - const DMReferenceDisplacementSample *displacement_sample = nullptr - ) const; - - [[nodiscard]] MappingStatus evaluate_face( - const ElementMappingData& element_data, - mfem::FaceElementTransformations& transformation, - FaceElementSide side, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - FaceMappingContext& context - ) const; - - [[nodiscard]] MappingStatus evaluate_point_variation( - const ElementMappingData& element_data, - const ElementDisplacementData& direction, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - const MappingPointContext& context, - DMWorkspace& workspace, - MappingPointVariation& variation - ) const; - - [[nodiscard]] MappingStatus evaluate_volume_variation( - const ElementMappingData& element_data, - const ElementDisplacementData& direction, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - const VolumeMappingContext& context, - DMWorkspace& workspace, - VolumeMappingVariation& variation - ) const; - - [[nodiscard]] MappingStatus evaluate_face_variation( - const ElementMappingData& element_data, - const ElementDisplacementData& direction, - mfem::FaceElementTransformations &transformation, - FaceElementSide side, - const mfem::IntegrationPoint& ip, - const FaceMappingContext& context, - DMWorkspace& workspace, - FaceMappingVariation& variation - ) const; - - [[nodiscard]] bool is_compactified_element(const mfem::ElementTransformation& transformation) const noexcept; - [[nodiscard]] int get_dimension() const noexcept; - [[nodiscard]] const compactification::ExteriorDomainMap& get_exterior_map() const noexcept; - private: - [[nodiscard]] MappingStatus evaluate_point_impl( - const ElementMappingData& element_data, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - MappingPointContext& context, - const DMFixedMeshGeometry* geometry, - const DMReferenceDisplacementSample* displacement_sample = nullptr - ) const; - - void validate_element_data( - const ElementMappingData& element_data - ) const; - - void evaluate_displacement_field( - const ElementDisplacementData& field, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - mfem::Vector& value, - mfem::DenseMatrix& jacobian, - const mfem::DenseMatrix *inverse_mesh_jacobian - ) const; - - [[nodiscard]] MappingStatus evaluate_compactification_coordinate( - const ElementCompactificationData& compactification, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - DMWorkspace& workspace, - CompactificationPointData& point_data, - const mfem::DenseMatrix *inverse_mesh_jacobian - ) const; - - [[nodiscard]] MappingStatus evaluate_point_variation_impl( - const ElementMappingData& element_data, - const ElementDisplacementData& direction, - mfem::ElementTransformation& transformation, - const mfem::IntegrationPoint& ip, - const MappingPointContext& context, - DMWorkspace& workspace, - MappingPointVariation& variation, - const mfem::DenseMatrix *inverse_mesh_jacobian - ) const; - - // These are static since they do not use any class members or workspace. Really these could and - // maybe should be free functions... - [[nodiscard]] static mfem::ElementTransformation& select_face_element_transformation( - const mfem::FaceElementTransformations& transformations, - FaceElementSide side - ); - - [[nodiscard]] static const mfem::IntegrationPoint& select_face_element_integration_point( - const mfem::FaceElementTransformations& transformations, - FaceElementSide side - ); - - private: - mutable std::unique_ptr m_geometryCache; - utils::user_options::DomainMapperOptions m_options; - std::unique_ptr m_exterior_map; - }; -} \ No newline at end of file diff --git a/src/include/serif/mapping/evaluator.hpp b/src/include/serif/mapping/evaluator.hpp deleted file mode 100644 index 80001d31..00000000 --- a/src/include/serif/mapping/evaluator.hpp +++ /dev/null @@ -1,93 +0,0 @@ -#pragma once -#include - -#include "mfem/fem/gridfunc.hpp" -#include "mfem/fem/eltrans.hpp" -#include "mfem/fem/intrules.hpp" -#include "mfem/linalg/vector.hpp" -#include "mfem/general/array.hpp" - -#include "serif/mapping/domain_mapper.hpp" -#include "serif/mapping/context.hpp" -#include "serif/mapping/status.hpp" - -namespace serif::mapping { - // Note this is NOT thread safe. - class GridFunctionMappingEvaluator { - public: - /* - * The evaluator references the supplied grid functions and caches copies of - * their element-local DOFs. You need to make sure you call InvalidateCache() or Refresh() after either - * grid function's values are modified. Changes to the finite-element-space sequence - * should be detected automatically. - * - * This object owns mutable workspace and cache state and *is not* thread-safe. - */ - GridFunctionMappingEvaluator( - const DomainMapper& mapper, - const mfem::GridFunction& displacement, - const mfem::GridFunction& compactification_coordinate - ); - - void invalidate_cache() noexcept; - void refresh(); - - [[nodiscard]] MappingStatus evaluate_point( - mfem::ElementTransformation &transformation, - const mfem::IntegrationPoint &integration_point, - MappingPointContext &context - ); - - [[nodiscard]] MappingStatus evaluate_volume( - mfem::ElementTransformation &transformation, - const mfem::IntegrationPoint &integration_point, - VolumeMappingContext &context - ); - - [[nodiscard]] MappingStatus evaluate_face( - mfem::FaceElementTransformations &transformation, - FaceElementSide side, - const mfem::IntegrationPoint &integration_point, - FaceMappingContext &context - ); - - [[nodiscard]] VolumeQuadratureContext get_quadrature_context( - mfem::ElementTransformation &transformation, - const mfem::IntegrationPoint &integration_point - ); - - [[nodiscard]] FaceQuadratureContext get_face_quadrature_context( - mfem::FaceElementTransformations &transformation, - const mfem::IntegrationPoint &integration_point, - FaceElementSide side = FaceElementSide::ELEMENT_1 - ); - - void get_physical_point( - mfem::ElementTransformation &transformation, - const mfem::IntegrationPoint &integration_point, - mfem::Vector &physical_position - ); - private: - void validate_field_bindings() const; - [[nodiscard]] bool invalidate_for_changed_spaces(); - void load_element(int element_id); - - private: - const DomainMapper &m_mapper; - const mfem::GridFunction &m_displacement; - const mfem::GridFunction &m_compactification_coordinate; - const mfem::FiniteElementSpace *m_displacement_space; - const mfem::FiniteElementSpace *m_compactification_space; - long m_displacement_space_sequence; - long m_compactification_space_sequence; - DMWorkspace m_workspace; - - mfem::Array m_displacement_dofs; - mfem::Array m_compactification_dofs; - mfem::Vector m_element_displacement; - mfem::Vector m_element_compactification; - std::unique_ptr m_displacement_data; - std::unique_ptr m_compactification_data; - int m_cached_element_id{-1}; - }; -} \ No newline at end of file diff --git a/src/include/serif/mapping/mapping.hpp b/src/include/serif/mapping/mapping.hpp new file mode 100644 index 00000000..0fda5fe8 --- /dev/null +++ b/src/include/serif/mapping/mapping.hpp @@ -0,0 +1,4 @@ +#pragma once + +#include "serif/backend/mfem/mapping/domain_mapper.hpp" +#include "serif/backend/mfem/mapping/transformations.hpp" diff --git a/src/include/serif/mapping/status.hpp b/src/include/serif/mapping/status.hpp deleted file mode 100644 index 3c08dc30..00000000 --- a/src/include/serif/mapping/status.hpp +++ /dev/null @@ -1,58 +0,0 @@ -#pragma once - -#include -#include -#include - -#include "serif/utils/error/codes.hpp" -#include "serif/meta/enums.hpp" - -namespace serif::mapping -{ - enum class MappingStatus : uint16_t { - VALID, - - // Start this at 400 so this is parallel to the mapping error codes. There must be a better way to - // do this but it is not worth the time to sort out now - INVALID_DIMENSION = 400, - NON_FINITE_MAPPING_INPUT, - INVALID_REFERENCE_RADIUS, - AT_COMPACTIFIED_INFINITY, - OUTSIDE_REFERENCE_DOMAIN, - NON_FINITE_MAPPING_RESULT, - NON_POSITIVE_DETERMINATE - }; - - static auto MappingStatus_enum_string_map = meta::enum_to_enum_string_map(); - static auto MappingStatus_string_enum_map = meta::enum_to_string_enum_map(); - - inline std::unordered_map make_MappingStatus_SERiFErrorCode_map() { - std::unordered_map ret; - for (const auto& [MS, MS_name] : MappingStatus_enum_string_map) { - if (std::to_underlying(MS) < 400) { - continue; // Skip VALID since it is not in the error code set - } - ret.insert({MS, utils::errors::SERiFErrorCodeNameToTypeMap.at(MS_name)}); - } - return ret; - } - - inline std::unordered_map make_SERIFErrorCode_MappingStatus_map() { - - std::unordered_map ret; - for (const auto& [MS, MS_name] : MappingStatus_enum_string_map) { - if (std::to_underlying(MS) < 400) { - continue; // Skip VALID since it is not in the error code set - } - ret.insert({utils::errors::SERiFErrorCodeNameToTypeMap.at(MS_name), MS}); - } - return ret; - } - - // We make these parallel lists so to ensure that we can lookup between SERiF error codes - // and equivalent MappingStatus. Note again this is ugly and bad and I am ashamed of this code.... - // anyhoodles there is likely a better way but pick your battles amiright? - static auto MappingStatus_SERiFErrorCode_map = make_MappingStatus_SERiFErrorCode_map(); - static auto SERiFErrorCode_MappingStatus_map = make_SERIFErrorCode_MappingStatus_map(); -} - diff --git a/src/include/serif/mapping/transformations.hpp b/src/include/serif/mapping/transformations.hpp deleted file mode 100644 index 02858ec6..00000000 --- a/src/include/serif/mapping/transformations.hpp +++ /dev/null @@ -1,87 +0,0 @@ -#pragma once - -#include "serif/mapping/variation.hpp" -#include "serif/mapping/context.hpp" - -namespace serif::mapping { - void map_hdiv_flux_to_physical( - const MappingPointContext& context, - const mfem::Vector& reference_flux, - mfem::Vector& physical_flux - ); - - void map_physical_flux_to_hdiv_reference( - const MappingPointContext& context, - const mfem::Vector& physical_flux, - mfem::Vector& reference_flux - ); - - void map_reference_gradient_to_physical( - const MappingPointContext& context, - const mfem::Vector& reference_gradient, - mfem::Vector& physical_gradient - ); - - void map_physical_gradient_to_reference( - const MappingPointContext& context, - const mfem::Vector& physical_gradient, - mfem::Vector& reference_gradient - ); - - void map_reference_vector_gradient_to_physical( - const MappingPointContext& context, - const mfem::DenseMatrix& reference_vector_gradient, - mfem::DenseMatrix& physical_vector_gradient - ); - - void map_physical_vector_gradient_to_reference( - const MappingPointContext& context, - const mfem::DenseMatrix& physical_vector_gradient, - mfem::DenseMatrix& reference_vector_gradient - ); - - [[nodiscard]] double map_hdiv_divergence_to_physical( - const MappingPointContext& context, - double reference_divergence - ); - - void compute_hdiv_mass_tensor( - const MappingPointContext& context, - mfem::DenseMatrix& mass_tensor - ); - - void compute_scalar_diffusion_tensor( - const MappingPointContext& context, - mfem::DenseMatrix& diffusion_tensor - ); - - void map_hcurl_field_to_physical( - const MappingPointContext& context, - const mfem::Vector& reference_field, - mfem::Vector& physical_field - ); - - void map_physical_field_to_hcurl_reference( - const MappingPointContext& context, - const mfem::Vector& physical_field, - mfem::Vector& reference_field - ); - - void map_hcurl_curl_to_physical( - const MappingPointContext& context, - const mfem::Vector &reference_curl, - mfem::Vector &physical_curl - ); - - void map_physical_curl_to_hcurl_reference( - const MappingPointContext& context, - const mfem::Vector& physical_curl, - mfem::Vector& reference_curl - ); - - void compute_hdiv_mass_tensor_variation( - const MappingPointContext& context, - const MappingPointVariation& variation, - mfem::DenseMatrix& mass_tensor_variation - ); -} \ No newline at end of file diff --git a/src/include/serif/mapping/variation.hpp b/src/include/serif/mapping/variation.hpp deleted file mode 100644 index 02fbdcd3..00000000 --- a/src/include/serif/mapping/variation.hpp +++ /dev/null @@ -1,29 +0,0 @@ -#pragma once - -#include "mfem/linalg/densemat.hpp" -#include "mfem/linalg/vector.hpp" - -namespace serif::mapping { - struct MappingPointVariation { - mfem::Vector displacement_variation; - mfem::Vector physical_position_variation; - mfem::DenseMatrix displacement_jacobian_variation; - mfem::DenseMatrix mapping_jacobian_variation; - mfem::DenseMatrix inverse_mapping_jacobian_variation; - double mapping_determinate_variation{0.0}; - }; - - struct VolumeMappingVariation { - MappingPointVariation mapping; - mfem::DenseMatrix inverse_element_jacobian_variation; - double weight_variation{0.0}; - }; - - struct FaceMappingVariation { - MappingPointVariation mapping; - mfem::Vector physical_normal_variation; - - double physical_surface_weight_variation{0.0}; - double normal_flux_scale_variation{0.0}; - }; -} \ No newline at end of file diff --git a/src/include/serif/mapping/volume_cache.hpp b/src/include/serif/mapping/volume_cache.hpp deleted file mode 100644 index 7c38105c..00000000 --- a/src/include/serif/mapping/volume_cache.hpp +++ /dev/null @@ -1,28 +0,0 @@ -#pragma once - -#include "mfem/linalg/densemat.hpp" - -#include "serif/mapping/context.hpp" -#include - -namespace serif::mapping { - class VolumeMappingCache { - public: - void set_size(int point_count, int dimension); - void store(int point, VolumeMappingContext& context); - void load(int point, VolumeMappingContext& context) const; - void load_inverse_jacobian(int point, mfem::DenseMatrix& inverse) const; - [[nodiscard]] int get_point_count() const noexcept; - [[nodiscard]] int get_dimension() const noexcept; - - private: - void validate_point(int point) const; - [[nodiscard]] const double* get_point_data(int point) const; - - std::vector m_data; - int m_point_count{0}; - int m_dimension{0}; - int m_point_stride{0}; - }; - -} \ No newline at end of file