From dc298a94c97622013dadd81be75ebe9eb4e4fa55 Mon Sep 17 00:00:00 2001 From: Nora Fayyazuddin Ljungberg Date: Wed, 1 May 2024 10:59:42 -0400 Subject: [PATCH] Correct tilde_b evolution in ValenciaDivClean --- docs/References.bib | 18 ++++++++++++++++++ .../Systems/GrMhd/ValenciaDivClean/Fluxes.cpp | 6 +++--- .../Systems/GrMhd/ValenciaDivClean/Fluxes.hpp | 2 +- .../Systems/GrMhd/ValenciaDivClean/Sources.cpp | 3 --- .../Systems/GrMhd/ValenciaDivClean/Sources.hpp | 2 +- .../Systems/GrMhd/ValenciaDivClean/System.hpp | 10 +++++++--- .../Systems/GrMhd/ValenciaDivClean/Fluxes.py | 2 +- .../Systems/GrMhd/ValenciaDivClean/Sources.py | 14 +++++--------- 8 files changed, 36 insertions(+), 21 deletions(-) diff --git a/docs/References.bib b/docs/References.bib index 0b569156b406..94b7bc61817f 100644 --- a/docs/References.bib +++ b/docs/References.bib @@ -2014,6 +2014,24 @@ @article{Palenzuela2010 year = {2010} } +@article{Palenzuela2018, + title = "A Simflowny-based finite-difference code for high-performance + computing in numerical relativity", + volume = "35", + ISSN = "1361-6382", + url = "http://dx.doi.org/10.1088/1361-6382/aad7f6", + DOI = "10.1088/1361-6382/aad7f6", + number = "18", + journal = "Classical and Quantum Gravity", + publisher = "IOP Publishing", + author = "Palenzuela, Carlos and Miñano, Borja and Viganò, Daniele and + Arbona, Antoni and Bona-Casas, Carles and Rigo, Andreu and + Bezares, Miguel and Bona, Carles and Massó, Joan", + year = "2018", + month = "aug", + pages="185007" +} + @article{Pareschi2005, author = {Pareschi, Lorenzo and Russo, Giovanni}, title = {Implicit-Explicit {Runge-Kutta} Schemes and Applications to diff --git a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.cpp b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.cpp index 2a4262f5b706..f2542dd44961 100644 --- a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.cpp +++ b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.cpp @@ -65,9 +65,9 @@ void fluxes_impl( tilde_s_flux->get(i, j) = tilde_s.get(j) * transport_velocity->get(i) - lapse_b_over_w.get(j) * tilde_b.get(i); tilde_b_flux->get(i, j) = - tilde_b.get(j) * transport_velocity->get(i) + - get(lapse) * (get(tilde_phi) * inv_spatial_metric.get(i, j) - - spatial_velocity.get(j) * tilde_b.get(i)); + tilde_b.get(j) * transport_velocity->get(i) - + transport_velocity->get(j) * tilde_b.get(i) + + get(lapse) * (get(tilde_phi) * inv_spatial_metric.get(i, j)); } tilde_s_flux->get(i, i) += get(pressure_star_lapse_sqrt_det_spatial_metric); } diff --git a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.hpp b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.hpp index 813a62585043..ddacddb53ff7 100644 --- a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.hpp +++ b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.hpp @@ -66,7 +66,7 @@ void fluxes_impl( * B^m v_m\\ * F^i({\tilde \tau}) = &~ {\tilde \tau} v^i_{tr} + \sqrt{\gamma} \alpha \left( * p + p_m \right) v^i - \alpha {\tilde B}^i B^m v_m \\ - * F^i({\tilde B}^j) = &~ {\tilde B}^j v^i_{tr} - \alpha v^j {\tilde B}^i + + * F^i({\tilde B}^j) = &~ {\tilde B}^j v^i_{tr} - v^j_{tr} {\tilde B}^i + * \alpha \gamma^{ij} {\tilde \Phi} \\ * F^i({\tilde \Phi}) = &~ \alpha {\tilde B^i} - \beta^i {\tilde \Phi} * \f} diff --git a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.cpp b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.cpp index a7fda42fa0ad..26e5a97b095f 100644 --- a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.cpp +++ b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.cpp @@ -157,9 +157,6 @@ void sources_impl( source_tilde_b->get(i) *= get(tilde_phi); source_tilde_b->get(i) -= get(lapse) * get(tilde_phi) * trace_spatial_christoffel_second.get(i); - for (size_t m = 0; m < 3; ++m) { - source_tilde_b->get(i) -= tilde_b.get(m) * d_shift.get(m, i); - } } trace(source_tilde_phi, extrinsic_curvature, inv_spatial_metric); diff --git a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.hpp b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.hpp index 212cd56e36a2..d8b9a4ce70b9 100644 --- a/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.hpp +++ b/src/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.hpp @@ -91,7 +91,7 @@ void sources_impl( * \alpha \\ * S({\tilde \tau}) = & \alpha {\tilde S}^{mn} K_{mn} - {\tilde S}^m \partial_m * \alpha \\ - * S({\tilde B}^i) = & -{\tilde B}^m \partial_m \beta^i + {\tilde \Phi} + * S({\tilde B}^i) = & {\tilde \Phi} * \gamma^{im} \partial_m \alpha + \alpha {\tilde \Phi} \left( \frac{1}{2} * \gamma^{il} \gamma^{jk} - \gamma^{ij} \gamma^{lk} \right) \partial_l * \gamma_{jk} \\ S({\tilde \Phi}) = & {\tilde B}^k \partial_k \alpha - \alpha K diff --git a/src/Evolution/Systems/GrMhd/ValenciaDivClean/System.hpp b/src/Evolution/Systems/GrMhd/ValenciaDivClean/System.hpp index 34c34874ddf4..3268b742726e 100644 --- a/src/Evolution/Systems/GrMhd/ValenciaDivClean/System.hpp +++ b/src/Evolution/Systems/GrMhd/ValenciaDivClean/System.hpp @@ -32,6 +32,8 @@ namespace grmhd { /// - Black hole-neutron star mergers with a hot nuclear equation of state: /// outflow and neutrino-cooled disk for a low-mass, high-spin case /// \cite Deaton2013 +/// - A Simflowny-based finite-difference code for high-performance computing +/// in Relativity \cite Palenzuela2018 /// namespace ValenciaDivClean { @@ -134,7 +136,7 @@ namespace ValenciaDivClean { * \delta^i_j - \alpha b_j \tilde{B}^i / W \\ * F^i({\tilde \tau}) &= {\tilde \tau} v^i_{tr} + \alpha \sqrt{\gamma} p^* v^i * - \alpha^2 b^0 \tilde{B}^i / W \\ - * F^i({\tilde B}^j) &= {\tilde B}^j v^i_{tr} - \alpha v^j {\tilde B}^i + * F^i({\tilde B}^j) &= {\tilde B}^j v^i_{tr} - v^j_{tr} {\tilde B}^i * + \alpha \gamma^{ij} {\tilde \Phi} \\ * F^i({\tilde \Phi}) &= \alpha {\tilde B^i} - {\tilde \Phi} \beta^i * \f} @@ -149,8 +151,7 @@ namespace ValenciaDivClean { * + {\tilde \tau}) \partial_j \alpha \\ * S({\tilde \tau}) &= \alpha {\tilde S}^{mn} K_{mn} * - {\tilde S}^m \partial_m \alpha \\ - * S({\tilde B}^j) &= -{\tilde B}^m \partial_m \beta^j - * + \Phi \partial_k (\alpha \sqrt{\gamma}\gamma^{jk}) \\ + * S({\tilde B}^j) &= \Phi \partial_k (\alpha \sqrt{\gamma}\gamma^{jk}) \\ * S({\tilde \Phi}) &= {\tilde B}^k \partial_k \alpha - \alpha K * {\tilde \Phi} - \alpha \kappa {\tilde \Phi} * \f} @@ -166,6 +167,9 @@ namespace ValenciaDivClean { * divergence-free (no-monopole) condition \f$\Phi = \partial_i {\tilde B}^i = * 0\f$. * + * \note The fluxes and sources of ${\tilde B}$ follow eq. 32 of + * \cite Palenzuela2018 rather than eq. 87 of \cite Moesta2014 + * * \note On the electron fraction side, the source term is currently set to * \f$S(\tilde{Y}_e) = 0\f$. Implementing the source term using neutrino scheme * is in progress (Last update : Oct 2022). diff --git a/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.py b/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.py index b89e6b560cba..af04f52697ec 100644 --- a/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.py +++ b/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Fluxes.py @@ -159,8 +159,8 @@ def tilde_b_flux( magnetic_field, ): result = np.outer(lapse * spatial_velocity - shift, tilde_b) + result -= np.outer(tilde_b, lapse * spatial_velocity - shift) result += lapse * inv_spatial_metric * tilde_phi - result -= lapse * np.outer(tilde_b, spatial_velocity) return result diff --git a/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.py b/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.py index c96cc0ace267..13823b0445e3 100644 --- a/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.py +++ b/tests/Unit/Evolution/Systems/GrMhd/ValenciaDivClean/Sources.py @@ -159,15 +159,11 @@ def source_tilde_b( ): term_one = np.einsum("ab, iab", inv_spatial_metric, d_spatial_metric) term_two = np.einsum("ab, aib", inv_spatial_metric, d_spatial_metric) - return ( - tilde_phi * np.einsum("a, ia", d_lapse, inv_spatial_metric) - - np.einsum("a, ai", tilde_b, d_shift) - + lapse - * tilde_phi - * ( - 0.5 * np.einsum("ia, a", inv_spatial_metric, term_one) - - np.einsum("ia, a", inv_spatial_metric, term_two) - ) + return tilde_phi * np.einsum( + "a, ia", d_lapse, inv_spatial_metric + ) + lapse * tilde_phi * ( + 0.5 * np.einsum("ia, a", inv_spatial_metric, term_one) + - np.einsum("ia, a", inv_spatial_metric, term_two) )