From 5ad65f26e0ac34eab2b8dbcef5ce80161ce852c5 Mon Sep 17 00:00:00 2001 From: hen-w Date: Tue, 22 Apr 2025 16:39:31 -0400 Subject: [PATCH] Generalize deriv_conformal_christoffel_second_kind() in CCZ4 Generalize the deriv_conformal_christoffel_second_kind() function within DerivChristoffel.hpp to support symmetric tensor parameter d_field_d for the second-order CCZ4 system --- src/DataStructures/Tensor/TypeAliases.hpp | 7 + .../Systems/Ccz4/DerivChristoffel.cpp | 130 +++++++++++------- .../Systems/Ccz4/DerivChristoffel.hpp | 16 ++- .../Systems/Ccz4/Test_DerivChristoffel.cpp | 10 ++ 4 files changed, 105 insertions(+), 58 deletions(-) diff --git a/src/DataStructures/Tensor/TypeAliases.hpp b/src/DataStructures/Tensor/TypeAliases.hpp index 9b226e9bd877..14835efc3d69 100644 --- a/src/DataStructures/Tensor/TypeAliases.hpp +++ b/src/DataStructures/Tensor/TypeAliases.hpp @@ -423,6 +423,13 @@ using iJkk = Tensor, SpatialIndex, SpatialIndex>>; +template +using iijj = Tensor, + index_list, + SpatialIndex, + SpatialIndex, + SpatialIndex>>; + template using iiJJ = Tensor, index_list, diff --git a/src/Evolution/Systems/Ccz4/DerivChristoffel.cpp b/src/Evolution/Systems/Ccz4/DerivChristoffel.cpp index bcdfc5bf0cf3..f3c0654d2de3 100644 --- a/src/Evolution/Systems/Ccz4/DerivChristoffel.cpp +++ b/src/Evolution/Systems/Ccz4/DerivChristoffel.cpp @@ -4,6 +4,7 @@ #include "Evolution/Systems/Ccz4/DerivChristoffel.hpp" #include +#include #include "DataStructures/DataVector.hpp" #include "DataStructures/Tensor/Tensor.hpp" @@ -11,46 +12,61 @@ #include "Utilities/Gsl.hpp" namespace Ccz4 { -template +template void deriv_conformal_christoffel_second_kind( const gsl::not_null*> result, const tnsr::II& inverse_conformal_spatial_metric, - const tnsr::ijj& field_d, - const tnsr::ijkk& d_field_d, + const tnsr::ijj& field_d, const TensorType& d_field_d, const tnsr::iJJ& field_d_up) { - for (size_t i = 0; i < Dim; ++i) { - for (size_t j = i; j < Dim; ++j) { - for (size_t k = 0; k < Dim; ++k) { - for (size_t m = 0; m < Dim; ++m) { - (*result).get(k, m, i, j) = - -2.0 * field_d_up.get(k, m, 0) * - (field_d.get(i, j, 0) + field_d.get(j, i, 0) - - field_d.get(0, i, j)) + - 0.5 * inverse_conformal_spatial_metric.get(m, 0) * - (d_field_d.get(k, i, j, 0) + d_field_d.get(i, k, j, 0) + - d_field_d.get(k, j, i, 0) + d_field_d.get(j, k, i, 0) - - d_field_d.get(k, 0, i, j) - d_field_d.get(0, k, i, j)); - for (size_t l = 1; l < Dim; ++l) { - (*result).get(k, m, i, j) += - -2.0 * field_d_up.get(k, m, l) * - (field_d.get(i, j, l) + field_d.get(j, i, l) - - field_d.get(l, i, j)) + - 0.5 * inverse_conformal_spatial_metric.get(m, l) * - (d_field_d.get(k, i, j, l) + d_field_d.get(i, k, j, l) + - d_field_d.get(k, j, i, l) + d_field_d.get(j, k, i, l) - - d_field_d.get(k, l, i, j) - d_field_d.get(l, k, i, j)); + // We keep this for loop specialization for faster speed when no symmetry + // of d_field_d is present. + if constexpr (std::is_same_v>) { + for (size_t i = 0; i < Dim; ++i) { + for (size_t j = i; j < Dim; ++j) { + for (size_t k = 0; k < Dim; ++k) { + for (size_t m = 0; m < Dim; ++m) { + (*result).get(k, m, i, j) = + -2.0 * field_d_up.get(k, m, 0) * + (field_d.get(i, j, 0) + field_d.get(j, i, 0) - + field_d.get(0, i, j)) + + 0.5 * inverse_conformal_spatial_metric.get(m, 0) * + (d_field_d.get(k, i, j, 0) + d_field_d.get(i, k, j, 0) + + d_field_d.get(k, j, i, 0) + d_field_d.get(j, k, i, 0) - + d_field_d.get(k, 0, i, j) - d_field_d.get(0, k, i, j)); + for (size_t l = 1; l < Dim; ++l) { + (*result).get(k, m, i, j) += + -2.0 * field_d_up.get(k, m, l) * + (field_d.get(i, j, l) + field_d.get(j, i, l) - + field_d.get(l, i, j)) + + 0.5 * inverse_conformal_spatial_metric.get(m, l) * + (d_field_d.get(k, i, j, l) + d_field_d.get(i, k, j, l) + + d_field_d.get(k, j, i, l) + d_field_d.get(j, k, i, l) - + d_field_d.get(k, l, i, j) - d_field_d.get(l, k, i, j)); + } } } } } + } else { + ::tenex::evaluate( + result, + -2.0 * field_d_up(ti::k, ti::M, ti::L) * + (field_d(ti::i, ti::j, ti::l) + field_d(ti::j, ti::i, ti::l) - + field_d(ti::l, ti::i, ti::j)) + + inverse_conformal_spatial_metric(ti::M, ti::L) * 0.5 * + (d_field_d(ti::k, ti::i, ti::j, ti::l) + + d_field_d(ti::i, ti::k, ti::j, ti::l) + + d_field_d(ti::k, ti::j, ti::i, ti::l) + + d_field_d(ti::j, ti::k, ti::i, ti::l) - + d_field_d(ti::k, ti::l, ti::i, ti::j) - + d_field_d(ti::l, ti::k, ti::i, ti::j))); } } -template +template tnsr::iJkk deriv_conformal_christoffel_second_kind( const tnsr::II& inverse_conformal_spatial_metric, - const tnsr::ijj& field_d, - const tnsr::ijkk& d_field_d, + const tnsr::ijj& field_d, const TensorType& d_field_d, const tnsr::iJJ& field_d_up) { tnsr::iJkk result{}; deriv_conformal_christoffel_second_kind(make_not_null(&result), @@ -94,47 +110,57 @@ deriv_contracted_conformal_christoffel_second_kind( #define DIM(data) BOOST_PP_TUPLE_ELEM(0, data) #define FRAME(data) BOOST_PP_TUPLE_ELEM(1, data) #define DTYPE(data) BOOST_PP_TUPLE_ELEM(2, data) +#define TTYPE(data) BOOST_PP_TUPLE_ELEM(3, data) -#define INSTANTIATE(_, data) \ +#define INSTANTIATE_CONFORMAL_CHRISTOFFEL(_, data) \ template void Ccz4::deriv_conformal_christoffel_second_kind( \ const gsl::not_null*> \ result, \ const tnsr::II& \ inverse_conformal_spatial_metric, \ const tnsr::ijj& field_d, \ - const tnsr::ijkk& d_field_d, \ + const TTYPE(data) < DTYPE(data), DIM(data), FRAME(data) > &d_field_d, \ const tnsr::iJJ& field_d_up); \ template tnsr::iJkk \ Ccz4::deriv_conformal_christoffel_second_kind( \ const tnsr::II& \ inverse_conformal_spatial_metric, \ const tnsr::ijj& field_d, \ - const tnsr::ijkk& d_field_d, \ - const tnsr::iJJ& field_d_up); \ - template void Ccz4::deriv_contracted_conformal_christoffel_second_kind( \ - const gsl::not_null*> \ - result, \ - const tnsr::II& \ - inverse_conformal_spatial_metric, \ - const tnsr::iJJ& field_d_up, \ - const tnsr::Ijj& \ - conformal_christoffel_second_kind, \ - const tnsr::iJkk& \ - d_conformal_christoffel_second_kind); \ - template tnsr::iJ \ - Ccz4::deriv_contracted_conformal_christoffel_second_kind( \ - const tnsr::II& \ - inverse_conformal_spatial_metric, \ - const tnsr::iJJ& field_d_up, \ - const tnsr::Ijj& \ - conformal_christoffel_second_kind, \ - const tnsr::iJkk& \ + const TTYPE(data) < DTYPE(data), DIM(data), FRAME(data) > &d_field_d, \ + const tnsr::iJJ& field_d_up); + +GENERATE_INSTANTIATIONS(INSTANTIATE_CONFORMAL_CHRISTOFFEL, (1, 2, 3), + (Frame::Grid, Frame::Inertial), (double, DataVector), + (tnsr::iijj, tnsr::ijkk)) + +#undef INSTANTIATE_CONFORMAL_CHRISTOFFEL +#undef TTYPE + +#define INSTANTIATE_CONTRACTED_CHRISTOFFEL(_, data) \ + template void Ccz4::deriv_contracted_conformal_christoffel_second_kind( \ + const gsl::not_null*> \ + result, \ + const tnsr::II& \ + inverse_conformal_spatial_metric, \ + const tnsr::iJJ& field_d_up, \ + const tnsr::Ijj& \ + conformal_christoffel_second_kind, \ + const tnsr::iJkk& \ + d_conformal_christoffel_second_kind); \ + template tnsr::iJ \ + Ccz4::deriv_contracted_conformal_christoffel_second_kind( \ + const tnsr::II& \ + inverse_conformal_spatial_metric, \ + const tnsr::iJJ& field_d_up, \ + const tnsr::Ijj& \ + conformal_christoffel_second_kind, \ + const tnsr::iJkk& \ d_conformal_christoffel_second_kind); -GENERATE_INSTANTIATIONS(INSTANTIATE, (1, 2, 3), (Frame::Grid, Frame::Inertial), - (double, DataVector)) +GENERATE_INSTANTIATIONS(INSTANTIATE_CONTRACTED_CHRISTOFFEL, (1, 2, 3), + (Frame::Grid, Frame::Inertial), (double, DataVector)) -#undef INSTANTIATE +#undef INSTANTIATE_CONTRACTED_CHRISTOFFEL #undef DTYPE #undef FRAME #undef DIM diff --git a/src/Evolution/Systems/Ccz4/DerivChristoffel.hpp b/src/Evolution/Systems/Ccz4/DerivChristoffel.hpp index 49efc27bd1a0..cd58f129558a 100644 --- a/src/Evolution/Systems/Ccz4/DerivChristoffel.hpp +++ b/src/Evolution/Systems/Ccz4/DerivChristoffel.hpp @@ -26,21 +26,25 @@ namespace Ccz4 { * `Ccz4::Tags::InverseConformalMetric`, the CCZ4 auxiliary variable defined by * `Ccz4::Tags::FieldD`, its spatial derivative, and the CCZ4 identity defined * by `Ccz4::Tags::FieldDUp`. + * \note In second-order Ccz4, we impose symmetry of index k and l + * in \f$ \partial_l D_{kij}=\frac{1}{2}\partial_l \partial_k + * \tilde{\gamma}_{ij} \f$, because partial derivatives commute and to use + * `second_partial_derivatives()`. \f$ D_{kij} \f$ is evolved in + * the first-order system so no such symmetry is imposed. */ -template +template void deriv_conformal_christoffel_second_kind( const gsl::not_null*> result, const tnsr::II& inverse_conformal_spatial_metric, - const tnsr::ijj& field_d, - const tnsr::ijkk& d_field_d, + const tnsr::ijj& field_d, const TensorType& d_field_d, const tnsr::iJJ& field_d_up); -template +template tnsr::iJkk deriv_conformal_christoffel_second_kind( const tnsr::II& inverse_conformal_spatial_metric, - const tnsr::ijj& field_d, - const tnsr::ijkk& d_field_d, + const tnsr::ijj& field_d, const TensorType& d_field_d, const tnsr::iJJ& field_d_up); + /// @} /// @{ diff --git a/tests/Unit/Evolution/Systems/Ccz4/Test_DerivChristoffel.cpp b/tests/Unit/Evolution/Systems/Ccz4/Test_DerivChristoffel.cpp index 6c12dbd1069e..1a383f3d103a 100644 --- a/tests/Unit/Evolution/Systems/Ccz4/Test_DerivChristoffel.cpp +++ b/tests/Unit/Evolution/Systems/Ccz4/Test_DerivChristoffel.cpp @@ -27,6 +27,16 @@ void test_compute_deriv_conformal_christoffel_second_kind( Frame::Inertial>), "DerivChristoffel", "deriv_conformal_christoffel_second_kind", {{{-1., 1.}}}, used_for_size); + pypp::check_with_random_values<1>( + static_cast (*)( + const tnsr::II&, + const tnsr::ijj&, + const tnsr::iijj&, + const tnsr::iJJ&)>( + &Ccz4::deriv_conformal_christoffel_second_kind), + "DerivChristoffel", "deriv_conformal_christoffel_second_kind", + {{{-1., 1.}}}, used_for_size); } template