Skip to content

Commit

Permalink
AsterSeeds: Add initial data for magnetic fields for BNS
Browse files Browse the repository at this point in the history
  • Loading branch information
jaykalinani committed Sep 3, 2024
1 parent b761341 commit a5e20b0
Show file tree
Hide file tree
Showing 5 changed files with 165 additions and 14 deletions.
9 changes: 6 additions & 3 deletions AsterSeeds/param.ccl
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@ KEYWORD test_case "Name of the testcase" STEERABLE=never

"spherical shock" :: ""
"magTOV" :: ""
"magBNS" :: ""
} "Balsara1"

private:
Expand Down Expand Up @@ -114,17 +115,19 @@ REAL r0 "Characteristic radial offset" STEERABLE=ALWAYS
*:* :: "Anything"
} 0.0

REAL dipole_x "x-coordinate of dipole center" STEERABLE=ALWAYS
# coordinates of dipole center, to be based on the location of maximum of rho

REAL dipole_x[2] "x-coordinate of dipole center" STEERABLE=ALWAYS
{
*:* :: "Anything"
} 0.0

REAL dipole_y "y-coordinate of the dipole center" STEERABLE=ALWAYS
REAL dipole_y[2] "y-coordinate of the dipole center" STEERABLE=ALWAYS
{
*:* :: "Anything"
} 0.0

REAL dipole_z "z-coordinate of the dipole center" STEERABLE=ALWAYS
REAL dipole_z[2] "z-coordinate of the dipole center" STEERABLE=ALWAYS
{
*:* :: "Anything"
} 0.0
Expand Down
26 changes: 24 additions & 2 deletions AsterSeeds/schedule.ccl
Original file line number Diff line number Diff line change
Expand Up @@ -53,19 +53,41 @@ if (CCTK_Equals(test_type, "3DTest")) {
if (CCTK_Equals(test_case, "magTOV")) {

STORAGE: Avec_cent
SCHEDULE AsterSeeds_InitializeCenteredAvec AT initial AFTER HydroBaseX_PostInitial BEFORE AsterX_InitialGroup {
SCHEDULE AsterSeeds_InitializeCenteredAvec_TOV AT initial AFTER HydroBaseX_PostInitial BEFORE AsterX_InitialGroup {
LANG: C
READS: HydroBaseX::press(everywhere)
WRITES: Avec_x_cent(everywhere) Avec_y_cent(everywhere) Avec_z_cent(everywhere)
}
"Set up initial conditions for the cell-centered vector potential"

SCHEDULE AsterSeeds_InitializeStagAvec AT initial AFTER AsterSeeds_InitializeCenteredAvec BEFORE AsterX_InitialGroup {
SCHEDULE AsterSeeds_InitializeStagAvec_TOV AT initial AFTER AsterSeeds_InitializeCenteredAvec_TOV BEFORE AsterX_InitialGroup {
LANG: C
READS: Avec_x_cent(everywhere) Avec_y_cent(everywhere) Avec_z_cent(everywhere)
WRITES: AsterX::Avec_x(interior) AsterX::Avec_y(interior) AsterX::Avec_z(interior)
SYNC: AsterX::Avec_x AsterX::Avec_y AsterX::Avec_z
}
"Set up initial conditions for the vector potential"
}

#Initial conditions for magnetic field for BNS

if (CCTK_Equals(test_case, "magBNS")) {

STORAGE: Avec_cent
SCHEDULE AsterSeeds_InitializeCenteredAvec_BNS AT initial AFTER HydroBaseX_PostInitial BEFORE AsterX_InitialGroup {
LANG: C
READS: HydroBaseX::press(everywhere)
WRITES: Avec_x_cent(everywhere) Avec_y_cent(everywhere) Avec_z_cent(everywhere)
}
"Set up initial conditions for the cell-centered vector potential"

SCHEDULE AsterSeeds_InitializeStagAvec_BNS AT initial AFTER AsterSeeds_InitializeCenteredAvec_BNS BEFORE AsterX_InitialGroup {
LANG: C
READS: Avec_x_cent(everywhere) Avec_y_cent(everywhere) Avec_z_cent(everywhere)
WRITES: AsterX::Avec_x(interior) AsterX::Avec_y(interior) AsterX::Avec_z(interior)
SYNC: AsterX::Avec_x AsterX::Avec_y AsterX::Avec_z
}
"Set up initial conditions for the vector potential"
}

}
118 changes: 118 additions & 0 deletions AsterSeeds/src/magBNS.cxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,118 @@
#include <loop_device.hxx>

#include <cctk.h>
#include <cctk_Arguments.h>
#include <cctk_Parameters.h>

#include <cmath>

#include "seeds_utils.hxx"

namespace AsterSeeds {
using namespace std;
using namespace Loop;
using namespace AsterUtils;

extern "C" void AsterSeeds_InitializeCenteredAvec_BNS(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeCenteredAvec_BNS;
DECLARE_CCTK_PARAMETERS;

if (CCTK_EQUALS(Afield_config, "internal dipole")) {

/* computing cell centered vector potential components */
grid.loop_all_device<1, 1, 1>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
// dummy initialization
CCTK_REAL x_local = 0.0;
CCTK_REAL y_local = 0.0;

// For star 1 at minus side
if (p.x < 0) {
x_local = p.x - dipole_x[0];
y_local = p.y - dipole_y[0];
}
// For star 2 at minus side
else {
x_local = p.x - dipole_x[1];
y_local = p.y - dipole_y[1];
}

CCTK_REAL Pcut = press_max * press_cut;
CCTK_REAL Pdiff = std::max(press(p.I) - Pcut, 0.0);
CCTK_REAL Aphi_local = Ab * pow(Pdiff, Avec_kappa);
Avec_x_cent(p.I) = -y_local * Aphi_local;
Avec_y_cent(p.I) = x_local * Aphi_local;
Avec_z_cent(p.I) = 0.0;
});

} else if (CCTK_EQUALS(Afield_config, "external dipole")) {

/* computing cell centered vector potential components */
grid.loop_all_device<1, 1, 1>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
// For star 1 at minus side
CCTK_REAL x_local_s1 = p.x - dipole_x[0];
CCTK_REAL y_local_s1 = p.y - dipole_y[0];
CCTK_REAL z_local_s1 = p.z - dipole_z[0];
CCTK_REAL cylrad2_s1 =
x_local_s1 * x_local_s1 + y_local_s1 * y_local_s1;
CCTK_REAL rsph_s1 =
sqrt(x_local_s1 * x_local_s1 + y_local_s1 * y_local_s1 +
z_local_s1 * z_local_s1);

CCTK_REAL Aphi_local_s1 =
B0 * (pow(r0, 3.0) / (pow(r0, 3.0) + pow(rsph_s1, 3.0))) /
sqrt(cylrad2_s1 + 1.0e-16);

// For star 2 at minus side
CCTK_REAL x_local_s2 = p.x - dipole_x[1];
CCTK_REAL y_local_s2 = p.y - dipole_y[1];
CCTK_REAL z_local_s2 = p.z - dipole_z[1];
CCTK_REAL cylrad2_s2 =
x_local_s2 * x_local_s2 + y_local_s2 * y_local_s2;
CCTK_REAL rsph_s2 =
sqrt(x_local_s2 * x_local_s2 + y_local_s2 * y_local_s2 +
z_local_s2 * z_local_s2);

CCTK_REAL Aphi_local_s2 =
B0 * (pow(r0, 3.0) / (pow(r0, 3.0) + pow(rsph_s2, 3.0))) /
sqrt(cylrad2_s2 + 1.0e-16);

Avec_x_cent(p.I) =
-(y_local_s1 * Aphi_local_s1 + y_local_s2 * Aphi_local_s2);
Avec_y_cent(p.I) =
x_local_s1 * Aphi_local_s1 + x_local_s2 * Aphi_local_s2;
Avec_z_cent(p.I) = 0.0;
});

} else {
CCTK_ERROR("Vector potential configuration not defined");
}
}

extern "C" void AsterSeeds_InitializeStagAvec_BNS(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeStagAvec_BNS;
DECLARE_CCTK_PARAMETERS;

grid.loop_int_device<1, 0, 0>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
Avec_x(p.I) = calc_avg_c2e(Avec_x_cent, p, 0);
});

grid.loop_int_device<0, 1, 0>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
Avec_y(p.I) = calc_avg_c2e(Avec_y_cent, p, 1);
});

grid.loop_int_device<0, 0, 1>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
Avec_z(p.I) = calc_avg_c2e(Avec_z_cent, p, 2);
});
}

} // namespace AsterSeeds
23 changes: 15 additions & 8 deletions AsterSeeds/src/magtov.cxx → AsterSeeds/src/magTOV.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -13,8 +13,8 @@ using namespace std;
using namespace Loop;
using namespace AsterUtils;

extern "C" void AsterSeeds_InitializeCenteredAvec(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeCenteredAvec;
extern "C" void AsterSeeds_InitializeCenteredAvec_TOV(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeCenteredAvec_TOV;
DECLARE_CCTK_PARAMETERS;

if (CCTK_EQUALS(Afield_config, "internal dipole")) {
Expand All @@ -23,11 +23,14 @@ extern "C" void AsterSeeds_InitializeCenteredAvec(CCTK_ARGUMENTS) {
grid.loop_all_device<1, 1, 1>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
CCTK_REAL x_local = p.x - dipole_x[0];
CCTK_REAL y_local = p.y - dipole_y[0];
CCTK_REAL z_local = p.z - dipole_z[0];
CCTK_REAL Pcut = press_max * press_cut;
CCTK_REAL Pdiff = std::max(press(p.I) - Pcut, 0.0);
CCTK_REAL Aphi_local = Ab * pow(Pdiff, Avec_kappa);
Avec_x_cent(p.I) = -p.y * Aphi_local;
Avec_y_cent(p.I) = p.x * Aphi_local;
Avec_x_cent(p.I) = -y_local * Aphi_local;
Avec_y_cent(p.I) = x_local * Aphi_local;
Avec_z_cent(p.I) = 0.0;
});

Expand All @@ -37,8 +40,12 @@ extern "C" void AsterSeeds_InitializeCenteredAvec(CCTK_ARGUMENTS) {
grid.loop_all_device<1, 1, 1>(
grid.nghostzones,
[=] CCTK_DEVICE(const Loop::PointDesc &p) CCTK_ATTRIBUTE_ALWAYS_INLINE {
CCTK_REAL cylrad2 = p.x * p.x + p.y * p.y;
CCTK_REAL rsph = sqrt(p.x * p.x + p.y * p.y + p.z * p.z);
CCTK_REAL x_local = p.x - dipole_x[0];
CCTK_REAL y_local = p.y - dipole_y[0];
CCTK_REAL z_local = p.z - dipole_z[0];
CCTK_REAL cylrad2 = x_local * x_local + y_local * y_local;
CCTK_REAL rsph =
sqrt(x_local * x_local + y_local * y_local + z_local * z_local);
CCTK_REAL rsph3 = pow(rsph, 3.0);
CCTK_REAL r03 = pow(r0, 3.0);
CCTK_REAL Aphi_local =
Expand All @@ -53,8 +60,8 @@ extern "C" void AsterSeeds_InitializeCenteredAvec(CCTK_ARGUMENTS) {
}
}

extern "C" void AsterSeeds_InitializeStagAvec(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeStagAvec;
extern "C" void AsterSeeds_InitializeStagAvec_TOV(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTSX_AsterSeeds_InitializeStagAvec_TOV;
DECLARE_CCTK_PARAMETERS;

grid.loop_int_device<1, 0, 0>(
Expand Down
3 changes: 2 additions & 1 deletion AsterSeeds/src/make.code.defn
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,8 @@ SRCS = \
2D_tests.cxx \
3D_tests.cxx \
atmosphere.cxx \
magtov.cxx
magTOV.cxx \
magBNS.cxx


# Subdirectories containing source files
Expand Down

0 comments on commit a5e20b0

Please sign in to comment.