Skip to content

Commit

Permalink
fixup
Browse files Browse the repository at this point in the history
  • Loading branch information
@iago-mendes
iago-mendes committed Aug 2, 2024
1 parent 20242c9 commit 72c568c
Show file tree
Hide file tree
Showing 4 changed files with 20 additions and 25 deletions.
22 changes: 10 additions & 12 deletions src/Elliptic/Systems/Xcts/Events/ObserveAdmIntegrals.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,6 @@ void local_adm_integrals(
const tnsr::II<DataVector, 3>& inv_spatial_metric,
const tnsr::ii<DataVector, 3>& extrinsic_curvature,
const Scalar<DataVector>& trace_extrinsic_curvature,
const tnsr::I<DataVector, 3, Frame::Inertial>& inertial_coords,
const InverseJacobian<DataVector, 3, Frame::ElementLogical,
Frame::Inertial>& inv_jacobian,
const Mesh<3>& mesh, const Element<3>& element,
Expand Down Expand Up @@ -81,8 +80,6 @@ void local_adm_integrals(
extrinsic_curvature, mesh, boundary_direction);
const auto face_trace_extrinsic_curvature = dg::project_tensor_to_boundary(
trace_extrinsic_curvature, mesh, boundary_direction);
const auto face_inertial_coords = dg::project_tensor_to_boundary(
inertial_coords, mesh, boundary_direction);
// This projection could be avoided by using
// domain::Tags::DetSurfaceJacobian from the DataBox, which is computed
// directly on the face (not projected). That would be better on Gauss
Expand All @@ -100,10 +97,16 @@ void local_adm_integrals(
face_inv_spatial_metric(ti::J, ti::L) *
face_extrinsic_curvature(ti::k, ti::l));

// Get interface mesh and normal
const auto& face_mesh = mesh.slice_away(boundary_direction.dimension());
const auto& conformal_face_normal =
conformal_face_normals.at(boundary_direction);

// Compute outward-pointing unit normal.
const auto face_r = magnitude(face_inertial_coords);
const auto face_unit_normal =
tenex::evaluate<ti::I>(face_inertial_coords(ti::I) / face_r());
// Note that `conformal_face_normal` is a one-form, so here we're simply
// transforming it into a vector.
const auto conformal_unit_normal = tenex::evaluate<ti::I>(
face_inv_conformal_metric(ti::I, ti::J) * conformal_face_normal(ti::j));

// Compute curved area element
const auto face_sqrt_det_conformal_metric =
Expand All @@ -112,11 +115,6 @@ void local_adm_integrals(
area_element(face_inv_jacobian, boundary_direction,
face_inv_conformal_metric, face_sqrt_det_conformal_metric);

// Get face mesh and face normal from data box
const auto& face_mesh = mesh.slice_away(boundary_direction.dimension());
const auto& conformal_face_normal =
conformal_face_normals.at(boundary_direction);

// Compute surface integrands
const auto mass_integrand = Xcts::adm_mass_surface_integrand(
face_deriv_conformal_factor, face_inv_conformal_metric,
Expand All @@ -128,7 +126,7 @@ void local_adm_integrals(
face_inv_extrinsic_curvature, face_trace_extrinsic_curvature);
const auto center_of_mass_integrand =
Xcts::center_of_mass_surface_integrand(face_conformal_factor,
face_unit_normal);
conformal_unit_normal);

// Contract surface integrands with face normal
const auto contracted_mass_integrand =
Expand Down
5 changes: 1 addition & 4 deletions src/Elliptic/Systems/Xcts/Events/ObserveAdmIntegrals.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,6 @@ void local_adm_integrals(
const tnsr::II<DataVector, 3>& inv_spatial_metric,
const tnsr::ii<DataVector, 3>& extrinsic_curvature,
const Scalar<DataVector>& trace_extrinsic_curvature,
const tnsr::I<DataVector, 3, Frame::Inertial>& inertial_coords,
const InverseJacobian<DataVector, 3, Frame::ElementLogical,
Frame::Inertial>& inv_jacobian,
const Mesh<3>& mesh, const Element<3>& element,
Expand Down Expand Up @@ -131,7 +130,6 @@ class ObserveAdmIntegrals : public Event {
gr::Tags::InverseSpatialMetric<DataVector, 3, Frame::Inertial>,
gr::Tags::ExtrinsicCurvature<DataVector, 3, Frame::Inertial>,
gr::Tags::TraceExtrinsicCurvature<DataVector>,
domain::Tags::Coordinates<3, Frame::Inertial>,
domain::Tags::InverseJacobian<3, Frame::ElementLogical, Frame::Inertial>,
domain::Tags::Mesh<3>, domain::Tags::Element<3>,
domain::Tags::Faces<3, domain::Tags::FaceNormal<3>>>;
Expand All @@ -149,7 +147,6 @@ class ObserveAdmIntegrals : public Event {
const tnsr::II<DataVector, 3>& inv_spatial_metric,
const tnsr::ii<DataVector, 3>& extrinsic_curvature,
const Scalar<DataVector>& trace_extrinsic_curvature,
const tnsr::I<DataVector, 3, Frame::Inertial>& inertial_coords,
const InverseJacobian<DataVector, 3, Frame::ElementLogical,
Frame::Inertial>& inv_jacobian,
const Mesh<3>& mesh, const Element<3>& element,
Expand All @@ -166,7 +163,7 @@ class ObserveAdmIntegrals : public Event {
deriv_conformal_factor, conformal_metric, inv_conformal_metric,
conformal_christoffel_second_kind, conformal_christoffel_contracted,
spatial_metric, inv_spatial_metric, extrinsic_curvature,
trace_extrinsic_curvature, inertial_coords, inv_jacobian, mesh, element,
trace_extrinsic_curvature, inv_jacobian, mesh, element,
conformal_face_normals);

// Save components of linear momentum as reduction data
Expand Down
4 changes: 2 additions & 2 deletions tests/Unit/Elliptic/Systems/Xcts/Events/Test_ObserveAdmIntegrals.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -253,8 +253,8 @@ void test_local_adm_integrals(const double& distance,
deriv_conformal_factor, conformal_metric, inv_conformal_metric,
conformal_christoffel_second_kind, conformal_christoffel_contracted,
spatial_metric, inv_spatial_metric, extrinsic_curvature,
trace_extrinsic_curvature, inertial_coords, inv_jacobian, mesh,
current_element, conformal_face_normals);
trace_extrinsic_curvature, inv_jacobian, mesh, current_element,
conformal_face_normals);
total_adm_mass.get() += get(local_adm_mass);
for (int I = 0; I < 3; I++) {
total_adm_linear_momentum.get(I) += local_adm_linear_momentum.get(I);
Expand Down
14 changes: 7 additions & 7 deletions tests/Unit/PointwiseFunctions/Xcts/Test_CenterOfMass.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@ void test_center_of_mass_surface_integral(const double distance) {
// Set up domain
const size_t h_refinement = 1;
const size_t p_refinement = 6;
domain::creators::Sphere shell{
const domain::creators::Sphere shell{
/* inner_radius */ z_shift + horizon_radius,
/* outer_radius */ distance,
/* interior */ domain::creators::Sphere::Excision{},
Expand Down Expand Up @@ -114,9 +114,9 @@ void test_center_of_mass_surface_integral(const double distance) {
shifted_fields);

// Compute outward-pointing unit normal.
const auto r = magnitude(inertial_coords);
const auto unit_normal =
tenex::evaluate<ti::I>(inertial_coords(ti::I) / r());
const auto conformal_r = magnitude(inertial_coords, conformal_metric);
const auto conformal_unit_normal =
tenex::evaluate<ti::I>(inertial_coords(ti::I) / conformal_r());

// Compute area element
const auto sqrt_det_conformal_metric =
Expand All @@ -126,8 +126,8 @@ void test_center_of_mass_surface_integral(const double distance) {
sqrt_det_conformal_metric);

// Integrate
const auto surface_integrand =
Xcts::center_of_mass_surface_integrand(conformal_factor, unit_normal);
const auto surface_integrand = Xcts::center_of_mass_surface_integrand(
conformal_factor, conformal_unit_normal);
for (int I = 0; I < 3; I++) {
surface_integral.get(I) += definite_integral(
surface_integrand.get(I) * get(conformal_area_element), face_mesh);
Expand All @@ -147,7 +147,7 @@ void test_center_of_mass_surface_integral(const double distance) {
SPECTRE_TEST_CASE("Unit.PointwiseFunctions.Xcts.CenterOfMass",
"[Unit][PointwiseFunctions]") {
// Test that integral converges with distance.
for (double distance : std::array<double, 3>({1.e3, 1.e4, 1.e5})) {
for (const double distance : std::array<double, 3>({1.e3, 1.e4, 1.e5})) {
test_center_of_mass_surface_integral(distance);
}
}

0 comments on commit 72c568c

Please sign in to comment.