From 2d7f5889b1d2197157aba0383d3f606ac0f805ce Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Thu, 5 Feb 2026 21:22:06 +0100 Subject: [PATCH 1/9] geometry: new InducedMetric class with g, g_inv, volume_element, and gradient --- lib/mito/geometry/induced_metric.h | 14 +++ .../geometry/induced_metric/InducedMetric.h | 100 ++++++++++++++++++ lib/mito/geometry/induced_metric/api.h | 19 ++++ lib/mito/geometry/induced_metric/factories.h | 25 +++++ lib/mito/geometry/induced_metric/forward.h | 28 +++++ lib/mito/geometry/induced_metric/public.h | 23 ++++ lib/mito/geometry/public.h | 3 + 7 files changed, 212 insertions(+) create mode 100644 lib/mito/geometry/induced_metric.h create mode 100644 lib/mito/geometry/induced_metric/InducedMetric.h create mode 100644 lib/mito/geometry/induced_metric/api.h create mode 100644 lib/mito/geometry/induced_metric/factories.h create mode 100644 lib/mito/geometry/induced_metric/forward.h create mode 100644 lib/mito/geometry/induced_metric/public.h diff --git a/lib/mito/geometry/induced_metric.h b/lib/mito/geometry/induced_metric.h new file mode 100644 index 00000000..efcf9556 --- /dev/null +++ b/lib/mito/geometry/induced_metric.h @@ -0,0 +1,14 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +// publish the induced metric interface +#include "induced/public.h" + + +// end of file diff --git a/lib/mito/geometry/induced_metric/InducedMetric.h b/lib/mito/geometry/induced_metric/InducedMetric.h new file mode 100644 index 00000000..6cced99b --- /dev/null +++ b/lib/mito/geometry/induced_metric/InducedMetric.h @@ -0,0 +1,100 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + template < + int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, + class ParametrizationT> + class InducedMetric { + + public: + // ambient dimension + static constexpr int D = coordsT::dim; + + // verify parametric coordinates have correct dimension + static_assert( + parametric_coordsT::dim == N, + "Parametric coordinates dimension must match submanifold dimension N"); + + // the metric space for the ambient coordinates + using ambient_metric_space_type = metric_space; + // the metric space for the parametric coordinates + using parametric_metric_space_type = metric_space; + // parametric coordinates type + using parametric_coordinates_type = parametric_coordsT; + // physical (ambient) coordinates + using coordinates_type = coordsT; + // induced metric tensor type (N x N) + using induced_metric_tensor_type = tensor::matrix_t; + // gradient vector type in ambient space + using gradient_vector_type = tensor::vector_t; + + private: + JacobianFieldT _jacobian; + ParametrizationT _parametrization; + + public: + constexpr InducedMetric(JacobianFieldT jacobian, ParametrizationT parametrization) : + _jacobian(std::move(jacobian)), + _parametrization(std::move(parametrization)) + {} + + public: + // Induced metric tensor at parametric point: g_induced = J^T * g * J + constexpr auto g(const parametric_coordinates_type & xi) const -> induced_metric_tensor_type + { + auto x = _parametrization(xi); + auto J = _jacobian(xi); + auto g_ambient = ambient_metric_space_type::g(x); + return tensor::transpose(J) * g_ambient * J; + } + + // Inverse induced metric at parametric point + constexpr auto g_inv(const parametric_coordinates_type & xi) const + -> induced_metric_tensor_type + { + return tensor::inverse(g(xi)); + } + + // Volume element: sqrt(det(g_induced)) + constexpr auto volume_element(const parametric_coordinates_type & xi) const + -> tensor::scalar_t + { + return std::sqrt(tensor::determinant(g(xi))); + } + + // Contravariant gradient in ambient space: J * g_induced^{-1} * dphi/dxi + // For N=1, dphi_dxi is a scalar; for N>1, dphi_dxi is vector_t + constexpr auto gradient( + const tensor::scalar_t & dphi_dxi, + const parametric_coordinates_type & xi) const -> gradient_vector_type + requires(N == 1) + { + auto J = _jacobian(xi); + auto g_inv_induced = g_inv(xi); + auto dphi_vec = tensor::vector_t<1>{ dphi_dxi }; + return J * (g_inv_induced * dphi_vec); + } + + constexpr auto gradient( + const tensor::vector_t & dphi_dxi, + const parametric_coordinates_type & xi) const -> gradient_vector_type + requires(N > 1) + { + auto J = _jacobian(xi); + auto g_inv_induced = g_inv(xi); + return J * (g_inv_induced * dphi_dxi); + } + }; + +} // namespace mito::geometry + + +// end of file diff --git a/lib/mito/geometry/induced_metric/api.h b/lib/mito/geometry/induced_metric/api.h new file mode 100644 index 00000000..cb44f8b7 --- /dev/null +++ b/lib/mito/geometry/induced_metric/api.h @@ -0,0 +1,19 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + // Factory: induced metric from Jacobian field and parametrization + template + constexpr auto induced_metric(JacobianT && jacobian, ParamT && parametrization); + +} // namespace mito::geometry + + +// end of file diff --git a/lib/mito/geometry/induced_metric/factories.h b/lib/mito/geometry/induced_metric/factories.h new file mode 100644 index 00000000..9b82ed12 --- /dev/null +++ b/lib/mito/geometry/induced_metric/factories.h @@ -0,0 +1,25 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + // Factory: create induced metric from Jacobian field and parametrization + template < + int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianT, + class ParamT> + constexpr auto induced_metric(JacobianT && jacobian, ParamT && parametrization) + { + return InducedMetric< + N, coordsT, parametric_coordsT, std::decay_t, std::decay_t>( + std::forward(jacobian), std::forward(parametrization)); + } +} // namespace mito::geometry + + +// end of file diff --git a/lib/mito/geometry/induced_metric/forward.h b/lib/mito/geometry/induced_metric/forward.h new file mode 100644 index 00000000..34cad59c --- /dev/null +++ b/lib/mito/geometry/induced_metric/forward.h @@ -0,0 +1,28 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + // class induced metric (N-dimensional submanifold in ambient coordsT) + template < + int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, + class ParametrizationT> + class InducedMetric; + + // alias for induced metric type + template < + int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, + class ParametrizationT> + using induced_metric_t = + InducedMetric; + +} + + +// end of file diff --git a/lib/mito/geometry/induced_metric/public.h b/lib/mito/geometry/induced_metric/public.h new file mode 100644 index 00000000..981d14a2 --- /dev/null +++ b/lib/mito/geometry/induced_metric/public.h @@ -0,0 +1,23 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +// get the forward declarations +#include "forward.h" + +// get the induced metric class +#include "InducedMetric.h" + +// get the factories +#include "factories.h" + +// get the api +#include "api.h" + + +// end of file diff --git a/lib/mito/geometry/public.h b/lib/mito/geometry/public.h index 9d67301b..f8878389 100644 --- a/lib/mito/geometry/public.h +++ b/lib/mito/geometry/public.h @@ -30,6 +30,9 @@ #include "polar.h" #include "spherical.h" +// induced metric (submanifolds) +#include "induced/public.h" + // classes implementation #include "CoordinateSystem.h" #include "Point.h" From 834dd879653059f3eb25d069fb7aa8398fdd6f36 Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Thu, 5 Feb 2026 21:25:51 +0100 Subject: [PATCH 2/9] geometry: add support for forms in the InducedMetric class --- .../geometry/induced_metric/InducedMetric.h | 56 +++++++++++++++++++ 1 file changed, 56 insertions(+) diff --git a/lib/mito/geometry/induced_metric/InducedMetric.h b/lib/mito/geometry/induced_metric/InducedMetric.h index 6cced99b..41cca199 100644 --- a/lib/mito/geometry/induced_metric/InducedMetric.h +++ b/lib/mito/geometry/induced_metric/InducedMetric.h @@ -92,6 +92,62 @@ namespace mito::geometry { auto g_inv_induced = g_inv(xi); return J * (g_inv_induced * dphi_dxi); } + + public: + // Basis elements for the parametric space (static, from parametric metric space) + template + requires(I >= 0 && I < N) + static constexpr auto e = parametric_metric_space_type::template e; + + private: + // Helper to create one-form from vector field and inverse induced metric + template + constexpr auto _dx() const + { + return functions::function([this](const parametric_coordinates_type & xi) { + // Get basis vector in parametric space + auto e_i = parametric_metric_space_type::template e(xi); + // Get inverse induced metric + auto g_inv_induced = this->g_inv(xi); + // Return the one-form: dx^i = g_inv * e_i + return tensor::one_form(e_i, g_inv_induced); + }); + } + + // Helper to wedge all N basis one-forms + template + constexpr auto _wedge(tensor::integer_sequence) const + requires(sizeof...(J) == N) + { + // Wedge all dx^J one-forms + return fields::wedge(_dx()...); + } + + public: + // The induced metric volume form (N-form on parametric domain) + // w = sqrt(det(g_induced)) * (dx^0 ∧ dx^1 ∧ ... ∧ dx^(N-1)) + constexpr auto w() const + { + // Create volume element field + auto volume_element_field = + functions::function([this](const parametric_coordinates_type & xi) { + return this->volume_element(xi); + }); + + // Wedge the basis one-forms + auto wedge_forms = _wedge(tensor::make_integer_sequence{}); + + // Return scaled wedge product + return volume_element_field * wedge_forms; + } + + // Basis one-forms (for convenience, if needed) + template + requires(I >= 0 && I < N) + constexpr auto dx() const + { + return _dx(); + } }; } // namespace mito::geometry From 7a03271855be59429513869cc13c7e06da871a29 Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Thu, 5 Feb 2026 21:36:12 +0100 Subject: [PATCH 3/9] geometry: fix includes after renaming the induced directory to induced_metric --- lib/mito/geometry/induced_metric.h | 2 +- lib/mito/geometry/public.h | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/lib/mito/geometry/induced_metric.h b/lib/mito/geometry/induced_metric.h index efcf9556..b672dfa2 100644 --- a/lib/mito/geometry/induced_metric.h +++ b/lib/mito/geometry/induced_metric.h @@ -8,7 +8,7 @@ // publish the induced metric interface -#include "induced/public.h" +#include "induced_metric/public.h" // end of file diff --git a/lib/mito/geometry/public.h b/lib/mito/geometry/public.h index f8878389..5554de0a 100644 --- a/lib/mito/geometry/public.h +++ b/lib/mito/geometry/public.h @@ -31,7 +31,7 @@ #include "spherical.h" // induced metric (submanifolds) -#include "induced/public.h" +#include "induced_metric/public.h" // classes implementation #include "CoordinateSystem.h" From c4f182eddba04565ee8165ede8e09b8232943d0f Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 6 Feb 2026 13:37:22 +0100 Subject: [PATCH 4/9] tests/geometry: tests for the induced metric calculations for segment and triangle and D from 1 to 3 --- .cmake/mito_tests_mito_lib.cmake | 2 + .../geometry/induced_metric_segment.cc | 167 ++++++++++++++++++ .../geometry/induced_metric_triangle.cc | 161 +++++++++++++++++ 3 files changed, 330 insertions(+) create mode 100644 tests/mito.lib/geometry/induced_metric_segment.cc create mode 100644 tests/mito.lib/geometry/induced_metric_triangle.cc diff --git a/.cmake/mito_tests_mito_lib.cmake b/.cmake/mito_tests_mito_lib.cmake index c4d6f232..5697eaba 100644 --- a/.cmake/mito_tests_mito_lib.cmake +++ b/.cmake/mito_tests_mito_lib.cmake @@ -31,6 +31,8 @@ mito_test_driver(tests/mito.lib/geometry/point.cc) mito_test_driver(tests/mito.lib/geometry/euclidean_metric_2D.cc) mito_test_driver(tests/mito.lib/geometry/euclidean_metric_3D.cc) mito_test_driver(tests/mito.lib/geometry/euclidean_submanifold_metric_3D.cc) +mito_test_driver(tests/mito.lib/geometry/induced_metric_segment.cc) +mito_test_driver(tests/mito.lib/geometry/induced_metric_triangle.cc) mito_test_driver(tests/mito.lib/geometry/cube_volume.cc) mito_test_driver(tests/mito.lib/geometry/metric.cc) mito_test_driver(tests/mito.lib/geometry/euclidean_metric_space.cc) diff --git a/tests/mito.lib/geometry/induced_metric_segment.cc b/tests/mito.lib/geometry/induced_metric_segment.cc new file mode 100644 index 00000000..e0412380 --- /dev/null +++ b/tests/mito.lib/geometry/induced_metric_segment.cc @@ -0,0 +1,167 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +#include +#include +#include + + +using coordinates_1d_t = mito::geometry::coordinates_t<1, mito::geometry::CARTESIAN>; +using coordinates_2d_t = mito::geometry::coordinates_t<2, mito::geometry::CARTESIAN>; +using coordinates_3d_t = mito::geometry::coordinates_t<3, mito::geometry::CARTESIAN>; + + +// helper function to compute directional derivative: grad · J_col +// this computes the inner product of the gradient with a tangent vector (column of Jacobian) +template +double +directional_derivative( + const mito::tensor::vector_t & grad, + const mito::tensor::matrix_t & J, int col_index) +{ + double result = 0.0; + for (int i = 0; i < AMBIENT_DIM; ++i) { + result += grad[i] * J[{ i, col_index }]; + } + return result; +} + +// helper function to compute gradient magnitude +template +double +gradient_magnitude(const mito::tensor::vector_t & grad) +{ + double result = 0.0; + for (int i = 0; i < DIM; ++i) { + result += grad[i] * grad[i]; + } + return std::sqrt(result); +} + + +TEST(Geometry, InducedMetricSegment1D) +{ + // unit segment in 1D: x(xi) = xi, J = 1, g_induced = 1, + // volume_element = 1 (equals length of segment) + auto jacobian = [](const coordinates_1d_t &) { + return mito::tensor::matrix_t<1, 1>{ 1.0 }; + }; + auto parametrization = [](const coordinates_1d_t & xi) { + return mito::geometry::cartesian::coordinates({ xi[0] }); + }; + + auto metric = mito::geometry::induced_metric<1, coordinates_1d_t, coordinates_1d_t>( + jacobian, parametrization); + + auto xi = mito::geometry::cartesian::coordinates<1>({ 0.5 }); + auto g_val = metric.g(xi); + EXPECT_DOUBLE_EQ(1.0, mito::tensor::trace(g_val)); // check g + EXPECT_DOUBLE_EQ(1.0, metric.volume_element(xi)); // check segment length + + // check inverse metric + auto g_inv = metric.g_inv(xi); + EXPECT_DOUBLE_EQ(1.0, mito::tensor::trace(g_inv)); + + // gradient: dphi/dxi = 1, check ⟨grad, tangent⟩ = dphi/dxi + auto grad = metric.gradient(1.0, xi); + auto J = jacobian(xi); + EXPECT_DOUBLE_EQ(1.0, directional_derivative(grad, J, 0)); + + // gradient magnitude: |grad| = 1/|tangent| = 1 + EXPECT_DOUBLE_EQ(1.0, gradient_magnitude(grad)); + + // nonlinear test: phi(xi) = xi^2, dphi/dxi = 2*xi + auto grad_nonlinear = metric.gradient(2.0 * xi[0], xi); + EXPECT_DOUBLE_EQ(2.0 * xi[0], directional_derivative(grad_nonlinear, J, 0)); + + // trigonometric test: phi(xi) = sin(xi), dphi/dxi = cos(xi) + auto grad_trig = metric.gradient(std::cos(xi[0]), xi); + EXPECT_DOUBLE_EQ(std::cos(xi[0]), directional_derivative(grad_trig, J, 0)); +} + + +TEST(Geometry, InducedMetricSegment2D) +{ + // diagonal segment from (0,0) to (3,4): x(xi) = (3*xi, 4*xi), J = [3; 4], g_induced = 25, + // volume element = 5 (equals length of segment) + auto jacobian = [](const coordinates_1d_t &) { + return mito::tensor::matrix_t<2, 1>{ 3.0, 4.0 }; + }; + auto parametrization = [](const coordinates_1d_t & xi) { + return mito::geometry::cartesian::coordinates({ 3.0 * xi[0], 4.0 * xi[0] }); + }; + + auto metric = mito::geometry::induced_metric<1, coordinates_2d_t, coordinates_1d_t>( + jacobian, parametrization); + + auto xi = mito::geometry::cartesian::coordinates<1>({ 0.5 }); + auto g_val = metric.g(xi); + EXPECT_DOUBLE_EQ(25.0, mito::tensor::trace(g_val)); // check g + EXPECT_DOUBLE_EQ(5.0, metric.volume_element(xi)); // check segment length + + // check inverse metric: g * g_inv should be identity + auto g_inv = metric.g_inv(xi); + EXPECT_DOUBLE_EQ(1.0 / 25.0, mito::tensor::trace(g_inv)); + + // gradient: dphi/dxi = 1, check ⟨grad, tangent⟩ = dphi/dxi + auto grad = metric.gradient(1.0, xi); + auto J = jacobian(xi); + EXPECT_DOUBLE_EQ(1.0, directional_derivative(grad, J, 0)); + + // gradient magnitude: |grad| = 1/|tangent| = 1/5 + EXPECT_DOUBLE_EQ(0.2, gradient_magnitude(grad)); + + // nonlinear test: phi(xi) = xi^2, dphi/dxi = 2*xi + auto grad_nonlinear = metric.gradient(2.0 * xi[0], xi); + EXPECT_DOUBLE_EQ(2.0 * xi[0], directional_derivative(grad_nonlinear, J, 0)); + + // trigonometric test: phi(xi) = sin(xi), dphi/dxi = cos(xi) + auto grad_trig = metric.gradient(std::cos(xi[0]), xi); + EXPECT_DOUBLE_EQ(std::cos(xi[0]), directional_derivative(grad_trig, J, 0)); +} + + +TEST(Geometry, InducedMetricSegment3D) +{ + // segment in 3D from (0,0,0) to (1,1,1): J = [1;1;1], g_induced = 3, + // volume_element = sqrt(3) (equals length of segment) + auto jacobian = [](const coordinates_1d_t &) { + return mito::tensor::matrix_t<3, 1>{ 1.0, 1.0, 1.0 }; + }; + auto parametrization = [](const coordinates_1d_t & xi) { + return mito::geometry::cartesian::coordinates({ xi[0], xi[0], xi[0] }); + }; + + auto metric = mito::geometry::induced_metric<1, coordinates_3d_t, coordinates_1d_t>( + jacobian, parametrization); + + auto xi = mito::geometry::cartesian::coordinates<1>({ 0.25 }); + auto g_val = metric.g(xi); + EXPECT_DOUBLE_EQ(3.0, mito::tensor::trace(g_val)); // check g + EXPECT_DOUBLE_EQ(std::sqrt(3.0), metric.volume_element(xi)); // check segment length + + // check inverse metric + auto g_inv = metric.g_inv(xi); + EXPECT_DOUBLE_EQ(1.0 / 3.0, mito::tensor::trace(g_inv)); + + // gradient: dphi/dxi = 1, check ⟨grad, tangent⟩ = dphi/dxi + auto grad = metric.gradient(1.0, xi); + auto J = jacobian(xi); + EXPECT_DOUBLE_EQ(1.0, directional_derivative(grad, J, 0)); + + // gradient magnitude: |grad| = 1/|tangent| = 1/sqrt(3) + EXPECT_DOUBLE_EQ(1.0 / std::sqrt(3.0), gradient_magnitude(grad)); + + // nonlinear test: phi(xi) = xi^2, dphi/dxi = 2*xi + auto grad_nonlinear = metric.gradient(2.0 * xi[0], xi); + EXPECT_DOUBLE_EQ(2.0 * xi[0], directional_derivative(grad_nonlinear, J, 0)); + + // trigonometric test: phi(xi) = sin(xi), dphi/dxi = cos(xi) + auto grad_trig = metric.gradient(std::cos(xi[0]), xi); + EXPECT_DOUBLE_EQ(std::cos(xi[0]), directional_derivative(grad_trig, J, 0)); +} + + +// end of file diff --git a/tests/mito.lib/geometry/induced_metric_triangle.cc b/tests/mito.lib/geometry/induced_metric_triangle.cc new file mode 100644 index 00000000..a2c688bf --- /dev/null +++ b/tests/mito.lib/geometry/induced_metric_triangle.cc @@ -0,0 +1,161 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +#include +#include +#include + + +using coordinates_2d_t = mito::geometry::coordinates_t<2, mito::geometry::CARTESIAN>; +using coordinates_3d_t = mito::geometry::coordinates_t<3, mito::geometry::CARTESIAN>; +using parametric_2d_t = mito::geometry::coordinates_t<2, mito::geometry::CARTESIAN>; + + +// helper function to compute directional derivative: grad · J_col +// this computes the inner product of the gradient with a tangent vector (column of Jacobian) +template +double +directional_derivative( + const mito::tensor::vector_t & grad, + const mito::tensor::matrix_t & J, int col_index) +{ + double result = 0.0; + for (int i = 0; i < AMBIENT_DIM; ++i) { + result += grad[i] * J[{ i, col_index }]; + } + return result; +} + +// helper function to test gradient in all parametric directions +// checks that ⟨grad, tangent_i⟩ = dphi[i] for each parametric direction i +template +void +test_gradient_directions( + const mito::tensor::vector_t & dphi, + const mito::tensor::vector_t & grad, + const mito::tensor::matrix_t & J) +{ + for (int col = 0; col < PARAM_DIM; ++col) { + EXPECT_DOUBLE_EQ(dphi[col], directional_derivative(grad, J, col)); + } +} + + +TEST(Geometry, InducedMetricTriangle2D) +{ + // unit right triangle: vertices (0,0), (1,0), (0,1) + // parametrization: x = (xi1 - xi0, 1 - xi0), J = [[-1, 1], [-1, 0]] + // this gives g_induced = J^T * J = [[2, -1], [-1, 1]] + auto jacobian = [](const parametric_2d_t &) { + mito::tensor::matrix_t<2, 2> J; + J[{ 0, 0 }] = -1; + J[{ 0, 1 }] = 1; + J[{ 1, 0 }] = -1; + J[{ 1, 1 }] = 0; + return J; + }; + auto parametrization = [](const parametric_2d_t & xi) { + return mito::geometry::cartesian::coordinates({ xi[1] - xi[0], 1.0 - xi[0] }); + }; + + auto metric = mito::geometry::induced_metric<2, coordinates_2d_t, parametric_2d_t>( + jacobian, parametrization); + + auto xi = mito::geometry::cartesian::coordinates<2>({ 1.0 / 3.0, 1.0 / 3.0 }); + auto g_val = metric.g(xi); + EXPECT_DOUBLE_EQ(2.0, (g_val[{ 0, 0 }])); + EXPECT_DOUBLE_EQ(-1.0, (g_val[{ 0, 1 }])); + EXPECT_DOUBLE_EQ(-1.0, (g_val[{ 1, 0 }])); + EXPECT_DOUBLE_EQ(1.0, (g_val[{ 1, 1 }])); + + // check symmetry of metric tensor + EXPECT_DOUBLE_EQ((g_val[{ 0, 1 }]), (g_val[{ 1, 0 }])); + EXPECT_DOUBLE_EQ(1.0, metric.volume_element(xi)); + + // check volume element equals sqrt(det(g)): det = 2*1 - (-1)*(-1) = 1 + EXPECT_DOUBLE_EQ(std::sqrt(mito::tensor::determinant(g_val)), metric.volume_element(xi)); + + // check inverse metric: g_inv should be [[1, 1], [1, 2]] + auto g_inv = metric.g_inv(xi); + EXPECT_DOUBLE_EQ(1.0, (g_inv[{ 0, 0 }])); + EXPECT_DOUBLE_EQ(1.0, (g_inv[{ 0, 1 }])); + EXPECT_DOUBLE_EQ(1.0, (g_inv[{ 1, 0 }])); + EXPECT_DOUBLE_EQ(2.0, (g_inv[{ 1, 1 }])); + + // gradient: dphi = (1, 0), check ⟨grad, tangent_i⟩ = dphi[i] + auto J = jacobian(xi); + auto dphi = mito::tensor::vector_t<2>{ 1.0, 0.0 }; + auto grad = metric.gradient(dphi, xi); + test_gradient_directions(dphi, grad, J); + + // gradient: dphi = (0, 1) + auto dphi2 = mito::tensor::vector_t<2>{ 0.0, 1.0 }; + auto grad2 = metric.gradient(dphi2, xi); + test_gradient_directions(dphi2, grad2, J); + + // gradient: dphi = (2, 3) + auto dphi3 = mito::tensor::vector_t<2>{ 2.0, 3.0 }; + auto grad3 = metric.gradient(dphi3, xi); + test_gradient_directions(dphi3, grad3, J); +} + + +TEST(Geometry, InducedMetricTriangle3D) +{ + // tilted triangle in 3D: (0,0,0), (1,0,1), (1,1,1), area = 0.5*sqrt(2) + // directors: (1,0,1), (1,1,1); cross = (-1,0,1), norm = sqrt(2), area = 0.5*sqrt(2) + auto jacobian = [](const parametric_2d_t &) { + // columns: (1,0,1), (1,1,1) + return mito::tensor::matrix_t<3, 2>{ 1.0, 1.0, 0.0, 1.0, 1.0, 1.0 }; + }; + auto parametrization = [](const parametric_2d_t & xi) { + auto x0 = mito::tensor::vector_t<3>{ 0.0, 0.0, 0.0 }; + auto x1 = mito::tensor::vector_t<3>{ 1.0, 0.0, 1.0 }; + auto x2 = mito::tensor::vector_t<3>{ 1.0, 1.0, 1.0 }; + auto x = (1.0 - xi[0] - xi[1]) * x0 + xi[0] * x1 + xi[1] * x2; + return mito::geometry::cartesian::coordinates({ x[0], x[1], x[2] }); + }; + + auto metric = mito::geometry::induced_metric<2, coordinates_3d_t, parametric_2d_t>( + jacobian, parametrization); + + auto xi = mito::geometry::cartesian::coordinates<2>({ 0.25, 0.25 }); + auto vol = metric.volume_element(xi); + EXPECT_DOUBLE_EQ(std::sqrt(2.0), vol); + + // check volume element equals sqrt(det(g)) + auto g_val_3d = metric.g(xi); + EXPECT_DOUBLE_EQ(std::sqrt(mito::tensor::determinant(g_val_3d)), metric.volume_element(xi)); + + // check symmetry of metric tensor + EXPECT_DOUBLE_EQ((g_val_3d[{ 0, 1 }]), (g_val_3d[{ 1, 0 }])); + + // gradient: dphi = (1, 0), check ⟨grad, tangent_i⟩ = dphi[i] + auto J = jacobian(xi); + auto dphi1 = mito::tensor::vector_t<2>{ 1.0, 0.0 }; + auto grad1 = metric.gradient(dphi1, xi); + test_gradient_directions(dphi1, grad1, J); + + // gradient: dphi = (0, 1) + auto dphi2 = mito::tensor::vector_t<2>{ 0.0, 1.0 }; + auto grad2 = metric.gradient(dphi2, xi); + test_gradient_directions(dphi2, grad2, J); + + // gradient: dphi = (2, 3) + auto dphi3 = mito::tensor::vector_t<2>{ 2.0, 3.0 }; + auto grad3 = metric.gradient(dphi3, xi); + test_gradient_directions(dphi3, grad3, J); + + // partition of unity: sum of basis function gradients should be zero + auto dphi0 = mito::tensor::vector_t<2>{ -1.0, -1.0 }; + auto grad0 = metric.gradient(dphi0, xi); + auto sum = grad0 + grad1 + grad2; + EXPECT_DOUBLE_EQ(0.0, sum[0]); + EXPECT_DOUBLE_EQ(0.0, sum[1]); + EXPECT_DOUBLE_EQ(0.0, sum[2]); +} + + +// end of file From 69b1a5aade00f0cb1a4b4603194f0bcb9487d4bd Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 6 Feb 2026 16:10:48 +0100 Subject: [PATCH 5/9] geometry: fix volume form construction --- .../geometry/induced_metric/InducedMetric.h | 30 +++++-------------- 1 file changed, 7 insertions(+), 23 deletions(-) diff --git a/lib/mito/geometry/induced_metric/InducedMetric.h b/lib/mito/geometry/induced_metric/InducedMetric.h index 41cca199..b13193e0 100644 --- a/lib/mito/geometry/induced_metric/InducedMetric.h +++ b/lib/mito/geometry/induced_metric/InducedMetric.h @@ -99,28 +99,20 @@ namespace mito::geometry { requires(I >= 0 && I < N) static constexpr auto e = parametric_metric_space_type::template e; - private: - // Helper to create one-form from vector field and inverse induced metric + public: + // Basis one-forms in parametric space (standard Euclidean dual basis) template - constexpr auto _dx() const - { - return functions::function([this](const parametric_coordinates_type & xi) { - // Get basis vector in parametric space - auto e_i = parametric_metric_space_type::template e(xi); - // Get inverse induced metric - auto g_inv_induced = this->g_inv(xi); - // Return the one-form: dx^i = g_inv * e_i - return tensor::one_form(e_i, g_inv_induced); - }); - } + requires(I >= 0 && I < N) + static constexpr auto dx = parametric_metric_space_type::template dx; + private: // Helper to wedge all N basis one-forms template - constexpr auto _wedge(tensor::integer_sequence) const + static constexpr auto _wedge(tensor::integer_sequence) requires(sizeof...(J) == N) { // Wedge all dx^J one-forms - return fields::wedge(_dx()...); + return fields::wedge(dx...); } public: @@ -140,14 +132,6 @@ namespace mito::geometry { // Return scaled wedge product return volume_element_field * wedge_forms; } - - // Basis one-forms (for convenience, if needed) - template - requires(I >= 0 && I < N) - constexpr auto dx() const - { - return _dx(); - } }; } // namespace mito::geometry From fc0a91d4b19bed1706530785538f6150e30b7ecd Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 6 Feb 2026 16:13:34 +0100 Subject: [PATCH 6/9] tests/geometry: add volume form checks to the InducedMetric test files --- .../geometry/induced_metric_segment.cc | 22 +++++++++++++ .../geometry/induced_metric_triangle.cc | 32 +++++++++++++++++++ 2 files changed, 54 insertions(+) diff --git a/tests/mito.lib/geometry/induced_metric_segment.cc b/tests/mito.lib/geometry/induced_metric_segment.cc index e0412380..84cfd2c0 100644 --- a/tests/mito.lib/geometry/induced_metric_segment.cc +++ b/tests/mito.lib/geometry/induced_metric_segment.cc @@ -120,6 +120,19 @@ TEST(Geometry, InducedMetricSegment2D) // trigonometric test: phi(xi) = sin(xi), dphi/dxi = cos(xi) auto grad_trig = metric.gradient(std::cos(xi[0]), xi); EXPECT_DOUBLE_EQ(std::cos(xi[0]), directional_derivative(grad_trig, J, 0)); + + // volume form test: w(tangent) = volume_element * tangent + auto w = metric.w(); + auto tangent = mito::tensor::vector_t<1>{ 1.0 }; + EXPECT_DOUBLE_EQ(5.0, w(xi)(tangent)); // volume_element = 5 + + // test scaling property + auto tangent2 = mito::tensor::vector_t<1>{ 2.0 }; + EXPECT_DOUBLE_EQ(10.0, w(xi)(tangent2)); + + // test with negative tangent + auto tangent_neg = mito::tensor::vector_t<1>{ -1.0 }; + EXPECT_DOUBLE_EQ(-5.0, w(xi)(tangent_neg)); } @@ -161,6 +174,15 @@ TEST(Geometry, InducedMetricSegment3D) // trigonometric test: phi(xi) = sin(xi), dphi/dxi = cos(xi) auto grad_trig = metric.gradient(std::cos(xi[0]), xi); EXPECT_DOUBLE_EQ(std::cos(xi[0]), directional_derivative(grad_trig, J, 0)); + + // volume form test: w(tangent) = volume_element * tangent + auto w = metric.w(); + auto tangent = mito::tensor::vector_t<1>{ 1.0 }; + EXPECT_DOUBLE_EQ(std::sqrt(3.0), w(xi)(tangent)); // volume_element = sqrt(3) + + // test scaling property + auto tangent2 = mito::tensor::vector_t<1>{ 0.5 }; + EXPECT_DOUBLE_EQ(0.5 * std::sqrt(3.0), w(xi)(tangent2)); } diff --git a/tests/mito.lib/geometry/induced_metric_triangle.cc b/tests/mito.lib/geometry/induced_metric_triangle.cc index a2c688bf..95497f03 100644 --- a/tests/mito.lib/geometry/induced_metric_triangle.cc +++ b/tests/mito.lib/geometry/induced_metric_triangle.cc @@ -99,6 +99,25 @@ TEST(Geometry, InducedMetricTriangle2D) auto dphi3 = mito::tensor::vector_t<2>{ 2.0, 3.0 }; auto grad3 = metric.gradient(dphi3, xi); test_gradient_directions(dphi3, grad3, J); + + // volume form test: w should give signed area of parallelogram spanned by two vectors + auto w = metric.w(); + auto v1 = mito::tensor::vector_t<2>{ 1.0, 0.0 }; + auto v2 = mito::tensor::vector_t<2>{ 0.0, 1.0 }; + EXPECT_DOUBLE_EQ(1.0, w(xi)(v1, v2)); // volume_element = 1 + + // test antisymmetry: w(v2, v1) = -w(v1, v2) + EXPECT_DOUBLE_EQ(-1.0, w(xi)(v2, v1)); + + // test with scaled vectors + auto v1_scaled = mito::tensor::vector_t<2>{ 2.0, 0.0 }; + auto v2_scaled = mito::tensor::vector_t<2>{ 0.0, 3.0 }; + EXPECT_DOUBLE_EQ(6.0, w(xi)(v1_scaled, v2_scaled)); // should scale as 2*3 = 6 + + // test with linearly dependent vectors (area should be zero) + auto v3 = mito::tensor::vector_t<2>{ 1.0, 0.0 }; + auto v4 = mito::tensor::vector_t<2>{ 2.0, 0.0 }; + EXPECT_DOUBLE_EQ(0.0, w(xi)(v3, v4)); } @@ -155,6 +174,19 @@ TEST(Geometry, InducedMetricTriangle3D) EXPECT_DOUBLE_EQ(0.0, sum[0]); EXPECT_DOUBLE_EQ(0.0, sum[1]); EXPECT_DOUBLE_EQ(0.0, sum[2]); + + // volume form test: w(director0, director1) gives the area element + auto w = metric.w(); + auto director0 = mito::tensor::vector_t<2>{ 1.0, 0.0 }; + auto director1 = mito::tensor::vector_t<2>{ 0.0, 1.0 }; + EXPECT_DOUBLE_EQ(std::sqrt(2.0), w(xi)(director0, director1)); + + // test antisymmetry + EXPECT_DOUBLE_EQ(-std::sqrt(2.0), w(xi)(director1, director0)); + + // test with linearly dependent vectors + auto v_parallel = mito::tensor::vector_t<2>{ 2.0, 0.0 }; + EXPECT_DOUBLE_EQ(0.0, w(xi)(director0, v_parallel)); } From e83e5a4ef23908286c90f6d8f8cdd0d26f394f5f Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 13 Feb 2026 14:32:48 +0100 Subject: [PATCH 7/9] geometry: replace the Induced Metric class with a static struct called pullback_metric --- .../geometry/induced_metric/InducedMetric.h | 140 ------------------ lib/mito/geometry/induced_metric/api.h | 19 --- lib/mito/geometry/induced_metric/factories.h | 25 ---- lib/mito/geometry/induced_metric/forward.h | 28 ---- .../{induced_metric.h => pullback_metric.h} | 4 +- lib/mito/geometry/pullback_metric/forward.h | 19 +++ lib/mito/geometry/pullback_metric/metric.h | 26 ++++ .../public.h | 10 +- 8 files changed, 49 insertions(+), 222 deletions(-) delete mode 100644 lib/mito/geometry/induced_metric/InducedMetric.h delete mode 100644 lib/mito/geometry/induced_metric/api.h delete mode 100644 lib/mito/geometry/induced_metric/factories.h delete mode 100644 lib/mito/geometry/induced_metric/forward.h rename lib/mito/geometry/{induced_metric.h => pullback_metric.h} (56%) create mode 100644 lib/mito/geometry/pullback_metric/forward.h create mode 100644 lib/mito/geometry/pullback_metric/metric.h rename lib/mito/geometry/{induced_metric => pullback_metric}/public.h (57%) diff --git a/lib/mito/geometry/induced_metric/InducedMetric.h b/lib/mito/geometry/induced_metric/InducedMetric.h deleted file mode 100644 index b13193e0..00000000 --- a/lib/mito/geometry/induced_metric/InducedMetric.h +++ /dev/null @@ -1,140 +0,0 @@ -// -*- c++ -*- -// -// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved -// - -// code guard -#pragma once - - -namespace mito::geometry { - - template < - int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, - class ParametrizationT> - class InducedMetric { - - public: - // ambient dimension - static constexpr int D = coordsT::dim; - - // verify parametric coordinates have correct dimension - static_assert( - parametric_coordsT::dim == N, - "Parametric coordinates dimension must match submanifold dimension N"); - - // the metric space for the ambient coordinates - using ambient_metric_space_type = metric_space; - // the metric space for the parametric coordinates - using parametric_metric_space_type = metric_space; - // parametric coordinates type - using parametric_coordinates_type = parametric_coordsT; - // physical (ambient) coordinates - using coordinates_type = coordsT; - // induced metric tensor type (N x N) - using induced_metric_tensor_type = tensor::matrix_t; - // gradient vector type in ambient space - using gradient_vector_type = tensor::vector_t; - - private: - JacobianFieldT _jacobian; - ParametrizationT _parametrization; - - public: - constexpr InducedMetric(JacobianFieldT jacobian, ParametrizationT parametrization) : - _jacobian(std::move(jacobian)), - _parametrization(std::move(parametrization)) - {} - - public: - // Induced metric tensor at parametric point: g_induced = J^T * g * J - constexpr auto g(const parametric_coordinates_type & xi) const -> induced_metric_tensor_type - { - auto x = _parametrization(xi); - auto J = _jacobian(xi); - auto g_ambient = ambient_metric_space_type::g(x); - return tensor::transpose(J) * g_ambient * J; - } - - // Inverse induced metric at parametric point - constexpr auto g_inv(const parametric_coordinates_type & xi) const - -> induced_metric_tensor_type - { - return tensor::inverse(g(xi)); - } - - // Volume element: sqrt(det(g_induced)) - constexpr auto volume_element(const parametric_coordinates_type & xi) const - -> tensor::scalar_t - { - return std::sqrt(tensor::determinant(g(xi))); - } - - // Contravariant gradient in ambient space: J * g_induced^{-1} * dphi/dxi - // For N=1, dphi_dxi is a scalar; for N>1, dphi_dxi is vector_t - constexpr auto gradient( - const tensor::scalar_t & dphi_dxi, - const parametric_coordinates_type & xi) const -> gradient_vector_type - requires(N == 1) - { - auto J = _jacobian(xi); - auto g_inv_induced = g_inv(xi); - auto dphi_vec = tensor::vector_t<1>{ dphi_dxi }; - return J * (g_inv_induced * dphi_vec); - } - - constexpr auto gradient( - const tensor::vector_t & dphi_dxi, - const parametric_coordinates_type & xi) const -> gradient_vector_type - requires(N > 1) - { - auto J = _jacobian(xi); - auto g_inv_induced = g_inv(xi); - return J * (g_inv_induced * dphi_dxi); - } - - public: - // Basis elements for the parametric space (static, from parametric metric space) - template - requires(I >= 0 && I < N) - static constexpr auto e = parametric_metric_space_type::template e; - - public: - // Basis one-forms in parametric space (standard Euclidean dual basis) - template - requires(I >= 0 && I < N) - static constexpr auto dx = parametric_metric_space_type::template dx; - - private: - // Helper to wedge all N basis one-forms - template - static constexpr auto _wedge(tensor::integer_sequence) - requires(sizeof...(J) == N) - { - // Wedge all dx^J one-forms - return fields::wedge(dx...); - } - - public: - // The induced metric volume form (N-form on parametric domain) - // w = sqrt(det(g_induced)) * (dx^0 ∧ dx^1 ∧ ... ∧ dx^(N-1)) - constexpr auto w() const - { - // Create volume element field - auto volume_element_field = - functions::function([this](const parametric_coordinates_type & xi) { - return this->volume_element(xi); - }); - - // Wedge the basis one-forms - auto wedge_forms = _wedge(tensor::make_integer_sequence{}); - - // Return scaled wedge product - return volume_element_field * wedge_forms; - } - }; - -} // namespace mito::geometry - - -// end of file diff --git a/lib/mito/geometry/induced_metric/api.h b/lib/mito/geometry/induced_metric/api.h deleted file mode 100644 index cb44f8b7..00000000 --- a/lib/mito/geometry/induced_metric/api.h +++ /dev/null @@ -1,19 +0,0 @@ -// -*- c++ -*- -// -// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved -// - -// code guard -#pragma once - - -namespace mito::geometry { - - // Factory: induced metric from Jacobian field and parametrization - template - constexpr auto induced_metric(JacobianT && jacobian, ParamT && parametrization); - -} // namespace mito::geometry - - -// end of file diff --git a/lib/mito/geometry/induced_metric/factories.h b/lib/mito/geometry/induced_metric/factories.h deleted file mode 100644 index 9b82ed12..00000000 --- a/lib/mito/geometry/induced_metric/factories.h +++ /dev/null @@ -1,25 +0,0 @@ -// -*- c++ -*- -// -// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved -// - -// code guard -#pragma once - - -namespace mito::geometry { - - // Factory: create induced metric from Jacobian field and parametrization - template < - int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianT, - class ParamT> - constexpr auto induced_metric(JacobianT && jacobian, ParamT && parametrization) - { - return InducedMetric< - N, coordsT, parametric_coordsT, std::decay_t, std::decay_t>( - std::forward(jacobian), std::forward(parametrization)); - } -} // namespace mito::geometry - - -// end of file diff --git a/lib/mito/geometry/induced_metric/forward.h b/lib/mito/geometry/induced_metric/forward.h deleted file mode 100644 index 34cad59c..00000000 --- a/lib/mito/geometry/induced_metric/forward.h +++ /dev/null @@ -1,28 +0,0 @@ -// -*- c++ -*- -// -// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved -// - -// code guard -#pragma once - - -namespace mito::geometry { - - // class induced metric (N-dimensional submanifold in ambient coordsT) - template < - int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, - class ParametrizationT> - class InducedMetric; - - // alias for induced metric type - template < - int N, coordinates_c coordsT, coordinates_c parametric_coordsT, class JacobianFieldT, - class ParametrizationT> - using induced_metric_t = - InducedMetric; - -} - - -// end of file diff --git a/lib/mito/geometry/induced_metric.h b/lib/mito/geometry/pullback_metric.h similarity index 56% rename from lib/mito/geometry/induced_metric.h rename to lib/mito/geometry/pullback_metric.h index b672dfa2..c9e14bda 100644 --- a/lib/mito/geometry/induced_metric.h +++ b/lib/mito/geometry/pullback_metric.h @@ -7,8 +7,8 @@ #pragma once -// publish the induced metric interface -#include "induced_metric/public.h" +// publish the pullback metric interface (replaces induced metric) +#include "pullback_metric/public.h" // end of file diff --git a/lib/mito/geometry/pullback_metric/forward.h b/lib/mito/geometry/pullback_metric/forward.h new file mode 100644 index 00000000..df44c9a3 --- /dev/null +++ b/lib/mito/geometry/pullback_metric/forward.h @@ -0,0 +1,19 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + // primary template: pullback of a metric along a parametrization + template + struct pullback_metric {}; + +} + + +// end of file diff --git a/lib/mito/geometry/pullback_metric/metric.h b/lib/mito/geometry/pullback_metric/metric.h new file mode 100644 index 00000000..87b22ba1 --- /dev/null +++ b/lib/mito/geometry/pullback_metric/metric.h @@ -0,0 +1,26 @@ +// -*- c++ -*- +// +// Copyright (c) 2020-2026, the MiTo Authors, all rights reserved +// + +// code guard +#pragma once + + +namespace mito::geometry { + + // specialization: pullback of ambient metric along a parametrization field + template + struct pullback_metric { + static constexpr auto field(const fields::field_c auto & parametrization) + { + auto g_ambient = ambient_metricT::field(); + auto jacobian = functions::derivative(parametrization); + return functions::transpose(jacobian) * g_ambient(parametrization) * jacobian; + } + }; + +} + + +// end of file diff --git a/lib/mito/geometry/induced_metric/public.h b/lib/mito/geometry/pullback_metric/public.h similarity index 57% rename from lib/mito/geometry/induced_metric/public.h rename to lib/mito/geometry/pullback_metric/public.h index 981d14a2..16d40e73 100644 --- a/lib/mito/geometry/induced_metric/public.h +++ b/lib/mito/geometry/pullback_metric/public.h @@ -10,14 +10,8 @@ // get the forward declarations #include "forward.h" -// get the induced metric class -#include "InducedMetric.h" - -// get the factories -#include "factories.h" - -// get the api -#include "api.h" +// get the pullback metric specialization +#include "metric.h" // end of file From e96d80b84eb656bf6579bfd6727817bc97e9e1b9 Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 13 Feb 2026 14:39:14 +0100 Subject: [PATCH 8/9] tensor: add a factory for building matrices from column vectors --- lib/mito/tensor/factories.h | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/lib/mito/tensor/factories.h b/lib/mito/tensor/factories.h index 008b86ea..447c96ae 100644 --- a/lib/mito/tensor/factories.h +++ b/lib/mito/tensor/factories.h @@ -47,6 +47,18 @@ namespace mito::tensor { return matrix * vector * v; }); } + // pack N column vectors of dimension D into a D x N matrix + template + constexpr auto matrix_from_columns(const std::array, N> & columns) -> matrix_t + { + matrix_t result; + for (int j = 0; j < N; ++j) { + for (int i = 0; i < D; ++i) { + result[{ i, j }] = columns[j][i]; + } + } + return result; + } } From 297c4dae07532cc9eedef9f7f25d210850420219 Mon Sep 17 00:00:00 2001 From: Pawel024 Date: Fri, 13 Feb 2026 14:42:21 +0100 Subject: [PATCH 9/9] manifolds: store a metric field instead of a volume form in manifolds and use this for integration with the Riemannian formula --- lib/mito/manifolds/Manifold.h | 34 +++++++-------- lib/mito/manifolds/api.h | 12 +++--- lib/mito/manifolds/factories.h | 75 ++++------------------------------ lib/mito/manifolds/forward.h | 6 +-- 4 files changed, 34 insertions(+), 93 deletions(-) diff --git a/lib/mito/manifolds/Manifold.h b/lib/mito/manifolds/Manifold.h index 55fb6bee..4b66470d 100644 --- a/lib/mito/manifolds/Manifold.h +++ b/lib/mito/manifolds/Manifold.h @@ -9,15 +9,15 @@ namespace mito::manifolds { - template + template requires(cellT::dim == coordsT::dim) class Manifold { private: // typedef for node using node_type = cellT::node_type; - // the volume form type - using volume_form_type = volumeFormT; + // the metric field type + using metric_field_type = metricFieldT; // the physical dimension of the manifold (that is that of the cell) static constexpr int D = cellT::dim; // the dimension of the manifold (that is that of the cell) @@ -38,10 +38,10 @@ namespace mito::manifolds { public: constexpr Manifold( const mesh_type & mesh, const coordinate_system_type & coordinate_system, - volume_form_type volume_form) : + metric_field_type metric_field) : _mesh(mesh), _coordinate_system(coordinate_system), - _volume_form(volume_form) + _metric_field(metric_field) {} // destructor @@ -123,19 +123,21 @@ namespace mito::manifolds { } private: - // computes the volume of a cell - template - constexpr auto _volume(const cell_type & cell, tensor::integer_sequence) const + // computes the volume of a cell via Riemannian formula (1/N!) * sqrt(det(J^T g J)) + template + constexpr auto _volume(const cell_type & cell, tensor::integer_sequence) const -> tensor::scalar_t - requires(sizeof...(J) == N) + requires(sizeof...(Is) == N) { // get the director edges of this cell and the point where they stem from auto [point, directors] = mito::geometry::directors(cell, _coordinate_system); - // compute the volume of a N-order simplicial cell as (1/N!) times the volume form - // contracted with the cell directors - auto volume = 1.0 / mito::tensor::factorial() * _volume_form(point)(directors[J]...); - // all done - return volume; + // evaluate the ambient metric at the base point + auto g = _metric_field(point); + // assemble the Jacobian from the directors (D x N matrix) + auto J = mito::tensor::matrix_from_columns(directors); + // compute the volume as (1/N!) * sqrt(det(J^T g J)) + return 1.0 / mito::tensor::factorial() + * std::sqrt(mito::tensor::determinant(mito::tensor::transpose(J) * g * J)); } private: @@ -143,8 +145,8 @@ namespace mito::manifolds { const mesh_type & _mesh; // the coordinate system const coordinate_system_type & _coordinate_system; - // the volume form - volume_form_type _volume_form; + // the metric field + metric_field_type _metric_field; }; } // namespace mito diff --git a/lib/mito/manifolds/api.h b/lib/mito/manifolds/api.h index a899f221..60cd7b6e 100644 --- a/lib/mito/manifolds/api.h +++ b/lib/mito/manifolds/api.h @@ -10,8 +10,8 @@ namespace mito::manifolds { // manifold alias - template - using manifold_t = Manifold; + template + using manifold_t = Manifold; // factory of manifolds from a mesh and a coordinate system template @@ -19,12 +19,12 @@ namespace mito::manifolds { const mesh::mesh_t & mesh, const geometry::coordinate_system_t & coordinate_system); - // factory submanifold from a mesh, a coordinate system and set of normal fields - template - constexpr auto submanifold( + // factory of manifolds with explicit metric field + template + constexpr auto manifold( const mesh::mesh_t & mesh, const geometry::coordinate_system_t & coordinate_system, - const mito::fields::vector_field_c auto & normal_field); + const metricFieldT & metric_field); } diff --git a/lib/mito/manifolds/factories.h b/lib/mito/manifolds/factories.h index efbd08f4..cea24ca2 100644 --- a/lib/mito/manifolds/factories.h +++ b/lib/mito/manifolds/factories.h @@ -9,85 +9,24 @@ namespace mito::manifolds { - // factory manifold - template + // factory manifold with explicit metric field + template constexpr auto manifold( const mesh::mesh_t & mesh, const geometry::coordinate_system_t & coordinate_system, - const volumeFormT & volume_form) + const metricFieldT & metric_field) { - return manifold_t(mesh, coordinate_system, volume_form); + return manifold_t(mesh, coordinate_system, metric_field); } - // factory of manifolds from a mesh and a coordinate system + // factory of manifolds from a mesh and a coordinate system (works for any dimension) template constexpr auto manifold( const mesh::mesh_t & mesh, const geometry::coordinate_system_t & coordinate_system) { - // the mesh type - using mesh_type = mesh::mesh_t; - - // assert that the manifold is of the highest dimension - static_assert(mesh_type::dim == mesh_type::order); - - // the metric space type - using metric_space_type = geometry::metric_space; - - // get the metric volume form - constexpr auto volume_form = metric_space_type::w; - - // return a new manifold - return manifold(mesh, coordinate_system, volume_form); - } - - // factory of submanifolds from a mesh, a coordinate system and set of normal fields - template - constexpr auto submanifold( - const mesh::mesh_t & mesh, - const geometry::coordinate_system_t & coordinate_system, - const fieldsT &... normal_field) - { - // the mesh type - using mesh_type = mesh::mesh_t; - - // assert that the number of normal fields provided corresponds to the dimension of the - // physical space minus the order of the manifold - static_assert(mesh_type::dim - mesh_type::order == sizeof...(fieldsT)); - - // the metric space type - using metric_space_type = geometry::metric_space; - - // get the metric volume form - constexpr auto w = metric_space_type::w; - - // strip namespace from the placeholder for forms contractions - using mito::tensor::_; - - // assert that the submanifold is either a line or a surface - // the general case would require the expansion of the {_} for {mesh_type::order} times in - // the contraction of the metric volume form - static_assert(mesh_type::order == 1 || mesh_type::order == 2); - - // if the submanifold is a line (either in 2D or 3D) then use one placeholder {_} - if constexpr (mesh_type::order == 1) { - // the restriction of the metric volume form to the submanifold - auto wS = functions::function( - [w, normal_field...](const coordsT & x) { return w(x)(normal_field(x)..., _); }); - - // return the manifold - return manifold(mesh, coordinate_system, wS); - } - - // if the submanifold is a surface (in 3D) then use two placeholders {_, _} - if constexpr (mesh_type::order == 2) { - // the restriction of the metric volume element to the submanifold - auto wS = functions::function( - [w, normal_field...](const coordsT & x) { return w(x)(normal_field(x)..., _, _); }); - - // return the manifold - return manifold(mesh, coordinate_system, wS); - } + constexpr auto g = geometry::metric::field(); + return manifold(mesh, coordinate_system, g); } } diff --git a/lib/mito/manifolds/forward.h b/lib/mito/manifolds/forward.h index e2811808..9dd3a1ef 100644 --- a/lib/mito/manifolds/forward.h +++ b/lib/mito/manifolds/forward.h @@ -10,7 +10,7 @@ namespace mito::manifolds { // class manifold - template + template requires(cellT::dim == coordsT::dim) class Manifold; @@ -18,8 +18,8 @@ namespace mito::manifolds { template concept manifold_c = requires(F c) { // require that F only binds to {Manifold} specializations - []( - const Manifold &) { + []( + const Manifold &) { }(c); }; }