Skip to content

GK passive advection & twistshift fixes #1084

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
JunoRavin merged 22 commits into from
Aug 5, 2026
Merged

GK passive advection & twistshift fixes #1084

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
22 commits
Select commit Hold shift + click to select a range
771c19a
Add a new GK equation object to do passive advection, with a velocity…
manauref Jul 3, 2026
1056915
Change field_id in rt_gk_passive so that the omega_H dt doesn't show …
manauref Jul 3, 2026
9aff987
Sync advection speeds so that parallel runs come out correctly. Add 2…
manauref Jul 4, 2026
b3601fd
I think the changes to rect_grid_find_cell in a recent PR might've br…
manauref Jul 6, 2026
fa5ca03
Add TS BC unit test corresponding to figure 14 of the paper. Seems to…
manauref Jul 6, 2026
965ca1a
Protect against integrating over zero-volume regions in some sub-cell…
manauref Jul 7, 2026
5af001f
Add passive advection CBC test. Works fine with constant q, but the q…
manauref Jul 7, 2026
ce9b41d
Add a TS BC unit test with the same setup as the CBC reg test. We may…
manauref Jul 7, 2026
9473d77
Add an assert in twistshift BC updater that errors out of the change …
manauref Jul 9, 2026
f4a9122
Change the way we look for shifted cells in periodic copies along y i…
manauref Jul 9, 2026
c289cd3
Fix y extents in CBC input files so the wedge number is indeed an int…
manauref Jul 9, 2026
ae49bfd
Fix kernel selection in device code, and add some missing casts.
manauref Jul 10, 2026
b6cac76
Add device CBC unit test. Add bc_twistshift_new method, instead of in…
manauref Jul 10, 2026
a9a5395
Merge branch 'gk_passive' of https://github.com/ammarhakim/gkylzero i…
manauref Jul 10, 2026
ce4c128
Change signature of bc_twistshift_new so that structs are passed by a…
manauref Jul 10, 2026
531b649
Merge branch 'main' into gk_passive
manauref Aug 4, 2026
74ec869
Merge branch 'main' into gk_passive
manauref Aug 4, 2026
4ebcbc8
Increase resolution in d3d_iwl and tcv_iwl tests so they run (they ca…
manauref Aug 4, 2026
5a6498d
Increase another tolerance to 1e-9 in TS BC unit test so it passes on…
manauref Aug 4, 2026
d17066b
Fix mem leak in GK app where an array was not released at program's end.
manauref Aug 5, 2026
520e8b7
Update perlmutter CPU machine files to remove all possibility of GPU …
manauref Aug 5, 2026
2e9eba0
Reduce Nmu in TCV iwl test to reduce runtime.
manauref Aug 5, 2026
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
1 change: 1 addition & 0 deletions core/zero/gkyl_eqn_type.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -186,4 +186,5 @@ enum gkyl_gk_collisionless_type {
GKYL_GK_COLLISIONLESS_EM_BPAR, // Electromagnetic with B_par fluctuations.
GKYL_GK_COLLISIONLESS_EM, // Electromagnetic with B_perp and B_par fluctuations.
GKYL_GK_COLLISIONLESS_NEUTRAL, // Neutral collisionless terms.
GKYL_GK_COLLISIONLESS_PASSIVE, // Passive advection (for dev/debug).
};
3 changes: 1 addition & 2 deletions core/zero/gkyl_rect_grid.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ struct gkyl_rect_grid* gkyl_rect_grid_new(int ndim,
* @param pick_lower If point is on a cell boundary, pick the lower cell if true (per direction).
* @param known_index Any known indices of where the point is (<0 if not known).
* @param cell_index Pointer to cell indices.
* Asserts: point lies within cell(s) specified by knownIdx (if specified).
* Asserts: point lies within cell(s) specified by known_index (if specified).
*/
GKYL_CU_DH
void gkyl_rect_grid_find_cell(const struct gkyl_rect_grid *grid, const double *point,
Expand Down Expand Up @@ -158,4 +158,3 @@ bool gkyl_rect_grid_read(struct gkyl_rect_grid *grid, FILE *fp);
* @param grid Grid object.
*/
void gkyl_rect_grid_release(struct gkyl_rect_grid *grid);

3 changes: 1 addition & 2 deletions core/zero/rect_grid.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,7 @@ gkyl_rect_grid_find_cell(const struct gkyl_rect_grid *grid, const double *point,
}
for (int i=0; i<search_num; i++){
for (int j=0; j<2; j++){
new_index[i] = GKYL_MAX2(GKYL_MIN2(mid_index[i] + plusminus[j],cells[i]),0);
new_index[i] = GKYL_MAX2(GKYL_MIN2(mid_index[i] + plusminus[j],cells[search_dim[i]]),1);
if (is_in_cell(grid, point, new_index, dim_trans, known_index))
low_high_index[j*nDim+search_dim[i]] = new_index[i];
}
Expand Down Expand Up @@ -212,4 +212,3 @@ gkyl_rect_grid_release(struct gkyl_rect_grid *grid)
{
gkyl_free(grid);
}

24 changes: 12 additions & 12 deletions gyrokinetic/apps/gk_field_2x3x.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -259,10 +259,10 @@ gk_field_2x3x_add_TS_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
.shear_dir = 0, // shift varies with x.
.edge = GKYL_LOWER_EDGE,
.cdim = app->cdim,
.bcdir_ext_update_r = app->global_par_ext,
.bcdir_ext_update_r = &app->global_par_ext,
.num_ghost = ghost, // one ghost per config direction
.basis = app->basis,
.grid = app->grid,
.basis = &app->basis,
.grid = &app->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -271,7 +271,7 @@ gk_field_2x3x_add_TS_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
T_LU_lo.shift_func = app->gk_geom->parallel_lower_bc_shift_func;
T_LU_lo.shift_func_ctx = app->gk_geom->parallel_lower_bc_shift_ctx;
}
f->bc_ts_lo = gkyl_bc_twistshift_new(&T_LU_lo);
f->bc_ts_lo = gkyl_bc_twistshift_inew(&T_LU_lo);

// TS BC updater for upper edge.
struct gkyl_bc_twistshift_inp T_UL_up = {
Expand All @@ -280,10 +280,10 @@ gk_field_2x3x_add_TS_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
.shear_dir = 0, // shift varies with x.
.edge = GKYL_UPPER_EDGE,
.cdim = app->cdim,
.bcdir_ext_update_r = app->global_par_ext,
.bcdir_ext_update_r = &app->global_par_ext,
.num_ghost = ghost, // one ghost per config direction
.basis = app->basis,
.grid = app->grid,
.basis = &app->basis,
.grid = &app->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -292,7 +292,7 @@ gk_field_2x3x_add_TS_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
T_UL_up.shift_func = app->gk_geom->parallel_upper_bc_shift_func;
T_UL_up.shift_func_ctx = app->gk_geom->parallel_upper_bc_shift_ctx;
}
f->bc_ts_up = gkyl_bc_twistshift_new(&T_UL_up);
f->bc_ts_up = gkyl_bc_twistshift_inew(&T_UL_up);

long buff_sz = app->global_lower_ghost[par_dir].volume;
f->bc_buffer = mkarr(app->use_gpu, app->basis.num_basis, buff_sz);
Expand Down Expand Up @@ -356,10 +356,10 @@ gk_field_2x3x_add_IWL_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
.shear_dir = 0, // shift varies with x.
.edge = GKYL_LOWER_EDGE,
.cdim = app->cdim,
.bcdir_ext_update_r = app->global_par_ext_core,
.bcdir_ext_update_r = &app->global_par_ext_core,
.num_ghost = ghost, // one ghost per config direction
.basis = app->basis,
.grid = app->grid,
.basis = &app->basis,
.grid = &app->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -368,7 +368,7 @@ gk_field_2x3x_add_IWL_updaters(struct gkyl_gyrokinetic_app *app, struct gk_field
T_LU_lo.shift_func = app->gk_geom->parallel_lower_bc_shift_func;
T_LU_lo.shift_func_ctx = app->gk_geom->parallel_lower_bc_shift_ctx;
}
f->bc_ts_lo = gkyl_bc_twistshift_new(&T_LU_lo);
f->bc_ts_lo = gkyl_bc_twistshift_inew(&T_LU_lo);

long buff_sz = GKYL_MAX2(app->global_lower_ghost_par_sol.volume, app->global_lower_ghost_par_core.volume);
f->bc_buffer = mkarr(app->use_gpu, app->basis.num_basis, buff_sz);
Expand Down
32 changes: 16 additions & 16 deletions gyrokinetic/apps/gk_species.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -886,10 +886,10 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
.shear_dir = 0, // shift varies with x.
.edge = GKYL_LOWER_EDGE,
.cdim = cdim,
.bcdir_ext_update_r = gks->local_par_ext,
.bcdir_ext_update_r = &gks->local_par_ext,
.num_ghost = ghost,
.basis = gks->basis,
.grid = gks->grid,
.basis = &gks->basis,
.grid = &gks->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -899,7 +899,7 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
tsinp.shift_func_ctx = app->gk_geom->parallel_lower_bc_shift_ctx;
}

gks->bc_ts_lo = gkyl_bc_twistshift_new(&tsinp);
gks->bc_ts_lo = gkyl_bc_twistshift_inew(&tsinp);

}
else if ( (gks->lower_bc[d].type == GKYL_BC_GK_SPECIES_COPY) ||
Expand Down Expand Up @@ -946,10 +946,10 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
.shear_dir = 0, // shift varies with x.
.edge = GKYL_UPPER_EDGE,
.cdim = cdim,
.bcdir_ext_update_r = gks->local_par_ext,
.bcdir_ext_update_r = &gks->local_par_ext,
.num_ghost = ghost,
.basis = gks->basis,
.grid = gks->grid,
.basis = &gks->basis,
.grid = &gks->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -959,7 +959,7 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
tsinp.shift_func_ctx = app->gk_geom->parallel_upper_bc_shift_ctx;
}

gks->bc_ts_up = gkyl_bc_twistshift_new(&tsinp);
gks->bc_ts_up = gkyl_bc_twistshift_inew(&tsinp);
}
else if ( (gks->upper_bc[d].type == GKYL_BC_GK_SPECIES_COPY) ||
(gks->upper_bc[d].type == GKYL_BC_GK_SPECIES_ABSORB) ||
Expand Down Expand Up @@ -999,10 +999,10 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
.shear_dir = 0, // shift varies with x.
.edge = GKYL_LOWER_EDGE,
.cdim = cdim,
.bcdir_ext_update_r = gks->local_par_ext_core,
.bcdir_ext_update_r = &gks->local_par_ext_core,
.num_ghost = ghost,
.basis = gks->basis,
.grid = gks->grid,
.basis = &gks->basis,
.grid = &gks->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -1011,18 +1011,18 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
tsinp_lo.shift_func = app->gk_geom->parallel_lower_bc_shift_func;
tsinp_lo.shift_func_ctx = app->gk_geom->parallel_lower_bc_shift_ctx;
}
gks->bc_ts_lo = gkyl_bc_twistshift_new(&tsinp_lo);
gks->bc_ts_lo = gkyl_bc_twistshift_inew(&tsinp_lo);

struct gkyl_bc_twistshift_inp tsinp_up = {
.bc_dir = par_dir,
.shift_dir = 1, // y shift.
.shear_dir = 0, // shift varies with x.
.edge = GKYL_UPPER_EDGE,
.cdim = cdim,
.bcdir_ext_update_r = gks->local_par_ext_core,
.bcdir_ext_update_r = &gks->local_par_ext_core,
.num_ghost = ghost,
.basis = gks->basis,
.grid = gks->grid,
.basis = &gks->basis,
.grid = &gks->grid,
.use_gpu = app->use_gpu,
};
if (app->gk_geom->geometry_id == GKYL_GEOMETRY_TOKAMAK)
Expand All @@ -1031,7 +1031,7 @@ gk_species_init_dynamic(struct gkyl_gk *gk_app_inp, struct gkyl_gyrokinetic_app
tsinp_up.shift_func = app->gk_geom->parallel_upper_bc_shift_func;
tsinp_up.shift_func_ctx = app->gk_geom->parallel_upper_bc_shift_ctx;
}
gks->bc_ts_up = gkyl_bc_twistshift_new(&tsinp_up);
gks->bc_ts_up = gkyl_bc_twistshift_inew(&tsinp_up);
}

// Set function pointers.
Expand Down
5 changes: 4 additions & 1 deletion gyrokinetic/apps/gk_species_bflux.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -688,8 +688,11 @@ gk_species_bflux_init(struct gkyl_gyrokinetic_app *app, void *species,
bflux->eqns = gkyl_malloc(bflux->num_eqns*sizeof(struct gkyl_dg_eqn *));

int eqc = 0;
if (gk_s->collisionless.collisionless_id)
if (gk_s->collisionless.collisionless_id == GKYL_GK_COLLISIONLESS_PASSIVE)
bflux->eqns[eqc++] = gkyl_dg_updater_gyrokinetic_passive_acquire_eqn(gk_s->collisionless.passive_slvr);
else if (gk_s->collisionless.collisionless_id)
bflux->eqns[eqc++] = gkyl_dg_updater_gyrokinetic_acquire_eqn(gk_s->collisionless.slvr);

if (gk_s->anom_diff.anom_diff_id)
bflux->eqns[eqc++] = gkyl_dg_updater_gk_anomalous_diffusion_acquire_eqn(gk_s->anom_diff.slvr);

Expand Down
Loading
Loading