Skip to content

Shock: add mirror_domain param #362

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
JamesMcClung merged 7 commits into from
Jun 11, 2025
Merged

Shock: add mirror_domain param #362

Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
a835196
psc_shock: double and mirror domain instead of reflecting
JamesMcClung Jun 5, 2025
4fbf549
psc_shock: -BGK fields
JamesMcClung Jun 5, 2025
47c7e78
psc_shock: reenable proper checks
JamesMcClung Jun 5, 2025
ea0c83a
psc_shock: fix mirrored ez
JamesMcClung Jun 11, 2025
b76c7e6
psc_shock: +mirrored_domain param
JamesMcClung Jun 11, 2025
8016022
psc_shock: use mirrored_domain param
JamesMcClung Jun 11, 2025
f77176a
psc_shock: rename to mirror_domain
JamesMcClung Jun 11, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
45 changes: 24 additions & 21 deletions src/psc_shock.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,6 @@ using PscConfig = PscConfig1vbecDouble<Dim>;

// ----------------------------------------------------------------------

using BgkMfields = PscConfig::Mfields;
using MfieldsState = PscConfig::MfieldsState;
using Mparticles = PscConfig::Mparticles;
using Balance = PscConfig::Balance;
Expand Down Expand Up @@ -73,6 +72,8 @@ double len_z;

int out_interval;

bool mirror_domain;

} // namespace

// ======================================================================
Expand Down Expand Up @@ -131,6 +132,8 @@ void setupParameters(int argc, char** argv)
int n_writes = parsedParams.getOrDefault<int>("n_writes", 100);
out_interval = psc_params.nmax / n_writes;

mirror_domain = parsedParams.getOrDefault<bool>("mirror_domain", false);

std::ifstream src(path_to_params, std::ios::binary);
std::ofstream dst("params_record.txt", std::ios::binary);
dst << src.rdbuf();
Expand All @@ -147,17 +150,19 @@ void setupParameters(int argc, char** argv)
Grid_t* setupGrid()
{
// FIXME add a check to catch mismatch between Dim and n grid points early
auto domain =
Grid_t::Domain{{nx, ny, nz}, // n grid points
{len_x, len_y, len_z}, // physical lengths
{0, 0, 0}, // location of lower corner
{n_patches_x, n_patches_y, n_patches_z}}; // n patches

auto bc =
psc::grid::BC{{BND_FLD_PERIODIC, BND_FLD_CONDUCTING_WALL, BND_FLD_PERIODIC},
{BND_FLD_PERIODIC, BND_FLD_CONDUCTING_WALL, BND_FLD_PERIODIC},
{BND_PRT_PERIODIC, BND_PRT_REFLECTING, BND_PRT_PERIODIC},
{BND_PRT_PERIODIC, BND_PRT_REFLECTING, BND_PRT_PERIODIC}};
int y_mult = mirror_domain ? 2 : 1;
auto domain = Grid_t::Domain{
{nx, y_mult * ny, nz}, // n grid points
{len_x, y_mult * len_y, len_z}, // physical lengths
{0, 0, 0}, // location of lower corner
{n_patches_x, y_mult * n_patches_y, n_patches_z}}; // n patches

auto BND_FLD_Y = mirror_domain ? BND_FLD_PERIODIC : BND_FLD_CONDUCTING_WALL;
auto BND_PRT_Y = mirror_domain ? BND_PRT_PERIODIC : BND_PRT_REFLECTING;
auto bc = psc::grid::BC{{BND_FLD_PERIODIC, BND_FLD_Y, BND_FLD_PERIODIC},
{BND_FLD_PERIODIC, BND_FLD_Y, BND_FLD_PERIODIC},
{BND_PRT_PERIODIC, BND_PRT_Y, BND_PRT_PERIODIC},
{BND_PRT_PERIODIC, BND_PRT_Y, BND_PRT_PERIODIC}};

auto kinds = Grid_t::Kinds(NR_KINDS);
kinds[KIND_ELECTRON] = {-1.0, electron_mass, "e"};
Expand Down Expand Up @@ -198,10 +203,11 @@ void initializeParticles(Balance& balance, Grid_t*& grid_ptr, Mparticles& mprts)
double temperature =
np.kind == KIND_ION ? ion_temperature : electron_temperature;
np.n = 1.0;
double vy_coef = mirror_domain && crd[1] > len_y ? -1 : 1;
np.p = setup_particles.createMaxwellian(
{np.kind,
np.n,
{v_upstream_x, v_upstream_y, v_upstream_z},
{v_upstream_x, vy_coef * v_upstream_y, v_upstream_z},
{temperature, temperature, temperature},
np.tag});
};
Expand All @@ -226,10 +232,11 @@ void fillGhosts(MF& mfld, int compBegin, int compEnd)
void initializeFields(MfieldsState& mflds)
{
setupFields(mflds, [&](int component, double coords[3]) {
double e_coef = mirror_domain && coords[1] > len_y ? -1 : 1;
switch (component) {
case EX: return e_x;
case EY: return e_y;
case EZ: return e_z;
case EX: return e_coef * e_x;
case EY: return e_coef * e_y;
case EZ: return e_coef * e_z;
case HX: return b_x;
case HY: return b_y;
case HZ: return b_z;
Expand Down Expand Up @@ -257,9 +264,6 @@ static void run(int argc, char** argv)
auto& grid = *grid_ptr;
MfieldsState mflds{grid};
Mparticles mprts{grid};
BgkMfields phi{grid, 1, mflds.ibn()};
BgkMfields gradPhi{grid, 3, mflds.ibn()};
BgkMfields divGradPhi{grid, 1, mflds.ibn()};

// ----------------------------------------------------------------------
// Set up various objects needed to run this case
Expand All @@ -279,8 +283,7 @@ static void run(int argc, char** argv)
// -- Checks
ChecksParams checks_params{};
checks_params.gauss.check_interval = out_interval;
// checks_params.gauss.dump_always = true;
checks_params.gauss.err_threshold = 1e-5;
checks_params.continuity.check_interval = out_interval;

Checks checks{grid, MPI_COMM_WORLD, checks_params};

Expand Down