Skip to content

Commit 40d2a21

Browse files
committed
Added ppm-plm recon
1 parent 568ce3e commit 40d2a21

4 files changed

Lines changed: 109 additions & 3 deletions

File tree

docs/input.md

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -46,11 +46,13 @@ Parameter: `reconstruction` (string)
4646
- `plm` : piecewise linear (second order)
4747
- `ppm` : piecewise parabolic (third order)
4848
- `mixed_plm_ppm` : `plm` for hydro variables and `ppm` for magnetic fields
49+
- `mixed_ppm_plm` : `ppm` for hydro variables and `plm` for magnetic fields
4950
- `limo3` : LimO3 (third order)
5051
- `weno3` : WENO3 (third order)
5152
- `wenoz` : WENO-Z (third order but more accurate than WENO3)
5253

53-
Note, `ppm` and `wenoz` need at least three ghost zones (`parthenon/mesh/num_ghost`).
54+
Note, `ppm`, `mixed_plm_ppm`, `mixed_ppm_plm`, and `wenoz` need at least three
55+
ghost zones (`parthenon/mesh/num_ghost`).
5456

5557
#### Floors
5658

src/hydro/hydro.cpp

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -21,6 +21,7 @@
2121
#include "../recon/dc_simple.hpp"
2222
#include "../recon/limo3_simple.hpp"
2323
#include "../recon/mixed_plm_ppm.hpp"
24+
#include "../recon/mixed_ppm_plm.hpp"
2425
#include "../recon/plm_simple.hpp"
2526
#include "../recon/ppm_simple.hpp"
2627
#include "../recon/weno3_simple.hpp"
@@ -339,6 +340,9 @@ std::shared_ptr<StateDescriptor> Initialize(ParameterInput *pin) {
339340
} else if (recon_str == "mixed_plm_ppm") {
340341
recon = Reconstruction::mixed_plm_ppm;
341342
recon_need_nghost = 3;
343+
} else if (recon_str == "mixed_ppm_plm") {
344+
recon = Reconstruction::mixed_ppm_plm;
345+
recon_need_nghost = 3;
342346
} else if (recon_str == "limo3") {
343347
recon = Reconstruction::limo3;
344348
recon_need_nghost = 2;
@@ -416,6 +420,8 @@ std::shared_ptr<StateDescriptor> Initialize(ParameterInput *pin) {
416420
add_flux_fun<Fluid::glmmhd, Reconstruction::ppm, RiemannSolver::hlle>(flux_functions);
417421
add_flux_fun<Fluid::glmmhd, Reconstruction::mixed_plm_ppm, RiemannSolver::hlle>(
418422
flux_functions);
423+
add_flux_fun<Fluid::glmmhd, Reconstruction::mixed_ppm_plm, RiemannSolver::hlle>(
424+
flux_functions);
419425
add_flux_fun<Fluid::glmmhd, Reconstruction::weno3, RiemannSolver::hlle>(flux_functions);
420426
add_flux_fun<Fluid::glmmhd, Reconstruction::limo3, RiemannSolver::hlle>(flux_functions);
421427
add_flux_fun<Fluid::glmmhd, Reconstruction::wenoz, RiemannSolver::hlle>(flux_functions);
@@ -424,6 +430,8 @@ std::shared_ptr<StateDescriptor> Initialize(ParameterInput *pin) {
424430
add_flux_fun<Fluid::glmmhd, Reconstruction::ppm, RiemannSolver::hlld>(flux_functions);
425431
add_flux_fun<Fluid::glmmhd, Reconstruction::mixed_plm_ppm, RiemannSolver::hlld>(
426432
flux_functions);
433+
add_flux_fun<Fluid::glmmhd, Reconstruction::mixed_ppm_plm, RiemannSolver::hlld>(
434+
flux_functions);
427435
add_flux_fun<Fluid::glmmhd, Reconstruction::weno3, RiemannSolver::hlld>(flux_functions);
428436
add_flux_fun<Fluid::glmmhd, Reconstruction::limo3, RiemannSolver::hlld>(flux_functions);
429437
add_flux_fun<Fluid::glmmhd, Reconstruction::wenoz, RiemannSolver::hlld>(flux_functions);

src/main.hpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -33,10 +33,10 @@ enum {
3333
enum { IV1 = 1, IV2 = 2, IV3 = 3, IPR = 4 };
3434

3535
enum class RiemannSolver { undefined, none, hlle, llf, hllc, hlld };
36-
enum class Reconstruction { undefined, dc, plm, ppm, mixed_plm_ppm, wenoz, weno3, limo3 };
36+
enum class Reconstruction { undefined, dc, plm, ppm, mixed_plm_ppm, mixed_ppm_plm, wenoz, weno3, limo3}
3737
enum class Integrator { undefined, rk1, rk2, vl2, rk3 };
3838
enum class Fluid { undefined, euler, glmmhd };
39-
enum class Cooling { none, tabular };
39+
enum class Cooling { none, tabular };  
4040
enum class Conduction { none, isotropic, anisotropic };
4141
enum class ConductionCoeff { none, fixed, spitzer };
4242
enum class Viscosity { none, isotropic };

src/recon/mixed_ppm_plm.hpp

Lines changed: 96 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,96 @@
1+
//========================================================================================
2+
// AthenaXXX astrophysical plasma code
3+
// Copyright(C) 2020 James M. Stone <jmstone@ias.edu> and the Athena code team
4+
// Licensed under the 3-clause BSD License (the "LICENSE")
5+
//========================================================================================
6+
#ifndef RECON_MIXED_PPM_PLM_HPP_
7+
#define RECON_MIXED_PPM_PLM_HPP_
8+
//! \file reconstruct_mixed_ppm_plm.hpp
9+
// \brief Mixed reconstruction: PPM for hydro variables (mass, mom1-3, energy),
10+
// PLM for magnetic fields (B1,B2,B3). Requires that PPM() and PLM() are
11+
// already defined elsewhere in the same compilation unit.
12+
// This version only works with uniform mesh spacing
13+
14+
#include <parthenon/parthenon.hpp>
15+
16+
#include "plm_simple.hpp"
17+
#include "ppm_simple.hpp"
18+
19+
using parthenon::ScratchPad2D;
20+
21+
//----------------------------------------------------------------------------------------
22+
//! \fn Reconstruct<Reconstruction::mixed_ppm_plm, int DIR>()
23+
// \brief Wrapper for mixed PPM/PLM reconstruction. Uses PPM for variable indices
24+
// 0-4 (density, momenta, energy) and PLM for indices 5-7 (Bcc1..Bcc3).
25+
// Matches the style and structure of PLM and PPM reconstruction files.
26+
27+
template <Reconstruction recon, int XNDIR>
28+
KOKKOS_INLINE_FUNCTION
29+
typename std::enable_if<recon == Reconstruction::mixed_ppm_plm, void>::type
30+
Reconstruct(parthenon::team_mbr_t const &member, const int k, const int j,
31+
const int il, const int iu, const parthenon::VariablePack<Real> &q,
32+
ScratchPad2D<Real> &ql, ScratchPad2D<Real> &qr) {
33+
const auto nvar = q.GetDim(4);
34+
35+
// Variable index mapping
36+
constexpr int ppm_lo = 0; // density
37+
constexpr int ppm_hi = 4; // energy
38+
constexpr int plm_lo = 5; // Bcc1
39+
constexpr int plm_hi = 7; // Bcc3
40+
41+
for (auto n = 0; n < nvar; ++n) {
42+
parthenon::par_for_inner(member, il, iu, [&](const int i) {
43+
//------------------------------------------------------------------------
44+
// PPM branch: density, momentum1-3, energy
45+
//------------------------------------------------------------------------
46+
if (n >= ppm_lo && n <= ppm_hi) {
47+
if constexpr (XNDIR == parthenon::X1DIR) {
48+
PPM(q(n, k, j, i - 2), q(n, k, j, i - 1), q(n, k, j, i), q(n, k, j, i + 1),
49+
q(n, k, j, i + 2), ql(n, i + 1), qr(n, i));
50+
} else if constexpr (XNDIR == parthenon::X2DIR) {
51+
PPM(q(n, k, j - 2, i), q(n, k, j - 1, i), q(n, k, j, i), q(n, k, j + 1, i),
52+
q(n, k, j + 2, i), ql(n, i), qr(n, i));
53+
} else if constexpr (XNDIR == parthenon::X3DIR) {
54+
PPM(q(n, k - 2, j, i), q(n, k - 1, j, i), q(n, k, j, i), q(n, k + 1, j, i),
55+
q(n, k + 2, j, i), ql(n, i), qr(n, i));
56+
} else {
57+
PARTHENON_FAIL("Unknown direction for mixed PPM branch.");
58+
}
59+
return;
60+
}
61+
62+
//------------------------------------------------------------------------
63+
// PLM branch: magnetic fields B1,B2,B3
64+
//------------------------------------------------------------------------
65+
if (n >= plm_lo && n <= plm_hi) {
66+
if constexpr (XNDIR == parthenon::X1DIR) {
67+
PLM(q(n, k, j, i - 1), q(n, k, j, i), q(n, k, j, i + 1), ql(n, i + 1),
68+
qr(n, i));
69+
} else if constexpr (XNDIR == parthenon::X2DIR) {
70+
PLM(q(n, k, j - 1, i), q(n, k, j, i), q(n, k, j + 1, i), ql(n, i), qr(n, i));
71+
} else if constexpr (XNDIR == parthenon::X3DIR) {
72+
PLM(q(n, k - 1, j, i), q(n, k, j, i), q(n, k + 1, j, i), ql(n, i), qr(n, i));
73+
} else {
74+
PARTHENON_FAIL("Unknown direction for mixed PLM branch.");
75+
}
76+
return;
77+
}
78+
79+
//------------------------------------------------------------------------
80+
// Default fallback: PPM
81+
//------------------------------------------------------------------------
82+
if constexpr (XNDIR == parthenon::X1DIR) {
83+
PPM(q(n, k, j, i - 2), q(n, k, j, i - 1), q(n, k, j, i), q(n, k, j, i + 1),
84+
q(n, k, j, i + 2), ql(n, i + 1), qr(n, i));
85+
} else if constexpr (XNDIR == parthenon::X2DIR) {
86+
PPM(q(n, k, j - 2, i), q(n, k, j - 1, i), q(n, k, j, i), q(n, k, j + 1, i),
87+
q(n, k, j + 2, i), ql(n, i), qr(n, i));
88+
} else if constexpr (XNDIR == parthenon::X3DIR) {
89+
PPM(q(n, k - 2, j, i), q(n, k - 1, j, i), q(n, k, j, i), q(n, k + 1, j, i),
90+
q(n, k + 2, j, i), ql(n, i), qr(n, i));
91+
}
92+
}); // par_for_inner
93+
} // for n
94+
}
95+
96+
#endif // RECON_MIXED_PPM_PLM_HPP_

0 commit comments

Comments
 (0)