Skip to content

Problems

Ben Prather edited this page Oct 3, 2024 · 7 revisions

Not all problems are necessarily documented here! Consult the kharma/prob/ directory and problem.cpp for up-to-date details. Not all magnetic field initializations are necessarily included either, see kharma/prob/seed_B.{cpp,hpp} for details on those. Most problems are highly configurable -- problem implementations are generally quite short, and are the best place to learn about how to set up the shock tube/customized Kelvin-Helmholtz/etc you might want.

Fishbone-Moncrief Torus

Equilibrium torus solution from Fishbone & Moncrief. Rephrased for convenience to be a function of the torus inner edge and the radius of pressure maximum, rather than angular momentum (the conversion uses the function lfish_calc in fm_torus.hpp). Name: torus, implemented in fm_torus.cpp, parameter block <torus>. Parameters:

  • rin: inner edge of torus. Default: 6.0
  • rmax: radius of torus pressure maximum. Default: 12.0
  • tilt: tilt the torus by this many degrees about the positive y-axis. The torus is still a midplane Fishbone-Moncrief solution -- it is built as if it was constructed in the midplane, then rotated after the fact to the new angle.

Perturbations

The <perturbation> block can be used to add random perturbations to any initial conditions KHARMA generates. It can be used with any problem, but is most often used in torus initializations to seed accretion. Note that random values are drawn independently per-zone, so the PSD of perturbations actually changes as resolution grows -- as zones get smaller, increasing power goes into higher spatial frequencies, which dissipate immediately when the simulation begins to evolve. However, the degree and distribution of perturbation does not seem to matter much: in SANE simulations, it only serves to seed the magnetorotational instability (MRI), and in MADs it is functionally useless as not even the MRI plays a primary role in promoting accretion. Parameters:

  • u_jitter: range over which to perturb the internal energy. The exact operation is $u \rightarrow u \cdot (1 + uniform(-\frac{u_j}{2}, \frac{u_j}{2}))$. Default: 0.0
  • jitter_above_rho: only apply the perturbation to zones above this density. Used to exclude floors when perturbing tori. Default: value of rho_min_geom
  • rng_seed: Random seed for repeatability. Default: 31337, obviously
  • serial: run perturbations serially on the CPU for better repeatability. Default: false

Torus-specific magnetic fields

Several of the options for seeding the magnetic field only make sense for tori (you absolutely can use them with other problems, they just wouldn't make as much sense). We outline them briefly. As parameters in the <b_field> block, we have additionally:

  • type: sane, mad, mad_quadrupole, r3s3, r5s5, vertical. sane and mad are EHT-standard definitions, see Wong et al.. mad_quadrupole adds a cosine factor to mad to obtain a quadrupole, whereas r3s3 removes the exponential decay but is otherwise identical to mad. r5s5 increases the rate of decay of the potential, concentrating more field toward the disk inner edge. vertical are old experiments toward generating MADs, YMMV.
  • beta_min: the desired ratio between the maximum gas pressure and maximum magnetic pressure in the starting disk. Note this is not equal to the minimum value of $\beta$ over the domain, but is instead the ratio of maxima -- the shorthand beta_min stays in use for historical reasons.
  • norm: whether to normalize the magnetic field to respect beta_min. Default on for tori, but off for all other problems.
  • A0: sets normalization for the vertical field, as a constant factor in front of the calculated vector potential (one can also normalize the vertical field with beta_min).
  • min_A: vector potential below which no magnetic field is seeded. This value is subtracted uniformly from the density-weighted initial vector potential, and the resulting potential applied only if it is positive (see Wong et al.). Default: 0.2

Bondi Problem

This initializes with the axisymmetric Bondi-Hoyle-Lyttleton solution. While this solution should be stationary if unmagnetized or radially magnetized (and indeed, this is a regression test) much of interest can be done by seeding initially symmetric inflow with a vertical or other magnetic field. Name: bondi, implemented in bondi.cpp, parameter block <bondi>. Parameters:

  • mdot: accretion rate, in code units -- effectively sets the simulation density scale.
  • rs: radius of the sonic point.
  • set_inner_bound: set the inner boundary to the Bondi solution (rather than an outflow boundary condition). Default: false
  • set_outer_bound: set the outer boundary to the Bondi solution at each step. Basically similar to Dirichlet boundaries, with different machinery. Default: true
  • ur_frac: set the radial velocity as this multiple of the equilibrium radial velocity. Default: 1.0
  • uphi: azimuthal velocity for all starting material, as a proportion of the Keplerian velocity. Default: 0
  • r_in_bondi: innermost radius to initialize with the analytic solution. Generally set outside the horizon for spinning black holes, since the Bondi solution is computed in Boyer-Lindquist coordinates.
  • fill_interior: Fill the region interior to r_in_bondi with the same density as r_in_bondi. Otherwise it will be filled with floor values.
  • diffinit: Use density falling as $r^{-1}$, more like observed converged states with magnetic field.

In addition to the parameters in <bondi>, setting boundaries/inner_x1=dirichlet or boundaries/outer_x1=dirichlet will activate the normal KHARMA Dirichlet/constant boundary conditions. There should be no practical difference from setting boundaries with the Bondi solution directly, but it may provide a small speed boost.

Magnetic Field Seeding

Because KHARMA supports both cell-centered and face-centered magnetic fields, problem setups generally do not set their own magnetic field values. Rather, they communicate the field parameters to SeedBField in seed_B.cpp, which takes care of setting the magnetic field at cell faces or centers, via the vector potential or directly.

Regardless of the problem you choose, you can also set these magnetic field parameters. This is useful for adding magnetization to existing problems: Kelvin-Helmholtz, explosions, or whatever you want to throw together. This is done by setting a number of options in the <b_field> block:

  • type: constant, wave, or shock_tube (in addition to the torus-specific magnetic fields). As described below.
  • B10, B20, B30, amp_B1, amp_B2, amp_B3, amp2_B1, amp2_B2, amp2_B3, k1, k2, k3, phase: used as described below.

The constant field sets $B^N =$ BN0.

The wave field sets $B^N =$ amp_BN $\cdot \cos($ k1 $\cdot X^1 +$ k2 $\cdot X^2 +$ k3 $\cdot X^3 +$ phase $) +$ amp2_BN $\cdot \sin($ k1 $\cdot X^1 +$ k2 $\cdot X^2 +$ k3 $\cdot X^3 +$ phase $)$.

The shock_tube field sets $B^N =$ amp_BN for $X^1 &lt;$ phase, $B^N =$ amp2_BN for $X^1 &gt;$ phase.

Resize & Restart

One of KHARMA's unique features is the ability to resize existing simulations, by reading in old simulation data, interpolating it onto a new finer (or coarser) grid, and cleaning up the resulting divergence in the magnetic field. This feature is handled the the resize_restart problem, with an example parameter file in pars/restarts/resize_restart.par, and an example usage in the test tests/resize/run.sh.

Basically, any KHARMA restart file (prob_name.out1.NNNNN.phdf) can be converted using pyharm into a barebones single-grid format formerly used by iharm3D, using the command pyharm convert --to_restart <filename>. The single-grid format is just for implementation convenience, it makes interpolation easier. Eventually, KHARMA will support arbitrary resizing from its own restart files, via the problem resize_restart_kharma, which is currently more limited.

The filename is then a parameter in a new problem parameter file, resize_restart.par, along with any parameters for the new simulation (such as the new grid size -- single grids only in 2024.9, but SMR/AMR grids are supported in dev), and parameters for the magnetic field "cleanup." When starting with the new parameter file, KHARMA will read the single-grid restart file, interpolate it onto the new grid, and then start a global solve of the elliptic equation $\nabla^2 p = \nabla \cdot B$. We know what the existing field divergence $\nabla \cdot B$ is, so this is just the Poisson equation. After solving for $p$, we know that $\nabla \cdot (B - \nabla p) = 0$, and that $\nabla p$ has no curl component, so we can take $B \rightarrow B - \nabla p$ to eliminate the divergence of the magnetic field without otherwise modifying it. (If you didn't follow that there's a full explanation in Tóth -- this is the projection scheme in his nomenclature, though in our case it is applied only once, with constrained transport thereafter).

In addition to the extra steps here regarding file formats, one should also be careful of which parameters are set during and after the magnetic field cleanup. For example, it is common to set parthenon/time/nlim=0 so that the cleanup step immediately outputs a usable restart file at the new resolution. However, this must be overridden again when using said restart file!

Tests

Most of the rest of KHARMA's many problem setups are tests (everything except the three problems above, as well as resize_restart_kharma and gizmo, which are only used for multi-zone runs). We run down the list with references here. The canonical problem setup for each test, along with any additional context or references, can be found in the test's parameter file (the parameters themselves, in pars/, not the testing scripts in tests/).

MHD/GRMHD

  • mhdmodes: linear modes of MHD, as implemented in every MHD code everywhere.
  • explosion: magnetized high-pressure explosion from Komissarov section 7.3.
  • kelvin_helmholtz: Setup for Kelvin-Helmholtz instability, in particular from Lecoanet et al.
  • orszag_tang: "vortex" problem producing a current sheet, Orszag & Tang
  • bz_monopole: Monopole magnetic field in vacuum (or rather, floor values). Good for checking against Blandford-Znajek solution, hence the name.

Electron heating

All tests for the electron temperature package are from Ressler et al..

  • hubble: Hubble flow with explicit heating for convergence test, section 4.1.1
  • noh: Noh shock, again testing heating fraction & stability, section 4.1.2
  • driven_turbulence: test for electron heating from turbulent dissipation, section 4.1.3

Extended MHD

All extended MHD tests are from Chandra et al..

  • anisotropic_conduction
  • conducting_atmosphere
  • emhdmodes
  • emhdshock