Skip to content

Commit

Permalink
Allow scaling of kappa in core rings
Browse files Browse the repository at this point in the history
Done to try and suppress mesh seam fluctuations. This works successfully but the fluctuations appear to be driven my ion momentum after all
  • Loading branch information
mikekryjak committed Aug 23, 2024
1 parent 784d621 commit 28f5825
Show file tree
Hide file tree
Showing 2 changed files with 24 additions and 0 deletions.
2 changes: 2 additions & 0 deletions include/evolve_pressure.hxx
Original file line number Diff line number Diff line change
Expand Up @@ -95,6 +95,8 @@ private:
bool low_p_diffuse_perp; ///< Add artificial cross-field diffusion at low electron pressure?

Field3D kappa_par; ///< Parallel heat conduction coefficient
BoutReal core_kappa_factor; ///< Scaling factor on first core ring conduction
int core_kappa_scale_rings; ///< Number of rings to scale conduction by core_kappa_factor

Field3D source, final_source; ///< External pressure source
Field3D Sp; ///< Total pressure source
Expand Down
22 changes: 22 additions & 0 deletions src/evolve_pressure.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -80,6 +80,15 @@ EvolvePressure::EvolvePressure(std::string name, Options& alloptions, Solver* so
.doc("Flux limiter factor. < 0 means no limit. Typical is 0.2 for electrons, 1 for ions.")
.withDefault(-1.0);

core_kappa_factor = options["core_kappa_factor"]
.doc("Multiply kappa by this factor in the first core ring. Helpful in avoiding mesh seam fluctuations at low rtol.")
.withDefault(1.0);

core_kappa_scale_rings = options["core_kappa_scale_rings"]
.doc("Number of rings to scale kappa by core_kappa_factor. Default is 1.")
.withDefault<int>(1);


p_div_v = options["p_div_v"]
.doc("Use p*Div(v) form? Default, false => v * Grad(p) form")
.withDefault<bool>(false);
Expand Down Expand Up @@ -343,6 +352,19 @@ void EvolvePressure::finally(const Options& state) {
}
}


// Loop over first core ring and multiply kappa by core_kappa_factor
for (int ix = mesh->firstX(); ix < mesh->firstX() + core_kappa_scale_rings + 1; ix++) {
for (int iy = 0; iy < mesh->LocalNy ; iy++) {
for (int iz = 0; iz < mesh->LocalNz; iz++) {

if (mesh->periodicY(ix)) { // select only periodic, i.e. core region
kappa_par(ix, iy, iz) *= core_kappa_factor;
}
}
}
}

// Note: Flux through boundary turned off, because sheath heat flux
// is calculated and removed separately
ddt(P) += (2. / 3) * FV::Div_par_K_Grad_par(kappa_par, T, false);
Expand Down

0 comments on commit 28f5825

Please sign in to comment.