You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
- optically thin cooling based on tabulated cooling tables with either Townsend 2009 exact integration or operator-split subcycling
25
+
- tracer particles
20
26
- static and adaptive mesh refinement
21
27
- problem generators for
22
28
- linear waves
@@ -79,7 +85,7 @@ The following examples are a few standard cases.
79
85
80
86
Most simple configuration (only CPU, no MPI).
81
87
The `Kokkos_ARCH_...` parameter should be adjusted to match the target machine where AthenaPK will be executed.
82
-
A full list of architecture keywords is available on the [Kokkos wiki](https://kokkos.github.io/kokkos-core-wiki/keywords.html#architecture-keywords).
88
+
A full list of architecture keywords is available on the [Kokkos wiki](https://kokkos.org/kokkos-core-wiki/get-started/configuration-guide.html#keywords-arch).
83
89
84
90
# configure with enabling Intel Broadwell or similar architecture (AVX2) instructions
#rkl2_max_dt_ratio = 100.0 # limits the ratio between the parabolic and hyperbolic timesteps (only used for RKL2 operator split integrator)
81
+
#cfl = 1.0 # Additional safety factor applied in the caluclation of the diffusive timestep (used in both unsplit and RKL2 integration schemes). Defaults to hyperbolic cfl.
82
+
83
+
conduction = anisotropic # none (disabled), or isotropic, or anisotropic
84
+
conduction_coeff = fixed # alternative: spitzer
85
+
thermal_diff_coeff_code = 0.01 # fixed coefficent in code units (code_length^2/code_time)
86
+
#spitzer_cond_in_erg_by_s_K_cm = 4.6e7 # spitzer coefficient in cgs units (requires definition of a unit system)
87
+
#conduction_sat_phi = 0.3 # fudge factor to account for uncertainties in saturated fluxes
88
+
89
+
90
+
viscosity = none # none (disabled) or isotropic
91
+
viscosity_coeff = fixed
92
+
mom_diff_coeff_code = 0.25 # fixed coefficent of the kinetmatic viscosity in code units (code_length^2/code_time)
93
+
94
+
resistivity = none # none (disabled) or ohmic
95
+
resistivity_coeff = fixed
96
+
ohm_diff_coeff_code = 0.25 # fixed coefficent of the magnetic (ohmic) diffusivity code units (code_length^2/code_time)
97
+
```
98
+
(An)isotropic thermal conduction (with fixed or Spitzer coefficient), and isotropic viscosity and
99
+
resistivity with fixed coefficient are currently implemented.
100
+
They can be integrated in an unsplit manner or operator split using a second-order accurate RKL2
101
+
supertimestepping algorithm.
102
+
More details are described in the following.
103
+
104
+
#### Integrators
105
+
106
+
Diffusive processes can be integrated in either an unsplit
107
+
fashion (`diffusion/integrator=unsplit`) or operator split using a second-order accurate
108
+
RKL2 super timestepping algorithm (`diffusion/integrator=rkl2`) following [^M+14].
109
+
110
+
In the unsplit case, the diffusive processes are included at the end of every stage in
111
+
the main integration loop and the global timestep is limited accordingly.
112
+
A separate CFL can be set for the diffusive processes via `diffusion/cfl=...`, which
113
+
defaults to the hyperbolic value if not set.
114
+
115
+
In the RKL2 case, the global timestep is not limited by the diffusive processes by default.
116
+
However, as reported by [^V+17] a large number of stages
117
+
($`s \approx \sqrt(\Delta t_{hyp}/\Delta t_{par}) \geq 20`$) in the supertimestepping
118
+
(in combination with anisotropic, limited diffusion) may lead to a loss in accuracy, which
119
+
is why the difference between hyperbolic and parabolic timesteps can be limited by
120
+
`diffusion/rkl2_max_dt_ratio=...` and a warning is shown if the ratio is above 400.
121
+
122
+
[^M+14]:
123
+
C. D. Meyer, D. S. Balsara, and T. D. Aslam, “A stabilized Runge–Kutta–Legendre method for explicit super-time-stepping of parabolic and mixed equations,” Journal of Computational Physics, vol. 257, pp. 594–626, 2014, doi: https://doi.org/10.1016/j.jcp.2013.08.021.
124
+
125
+
[^V+17]:
126
+
B. Vaidya, D. Prasad, A. Mignone, P. Sharma, and L. Rickler, “Scalable explicit implementation of anisotropic diffusion with Runge–Kutta–Legendre super-time stepping,” Monthly Notices of the Royal Astronomical Society, vol. 472, no. 3, pp. 3147–3160, 2017, doi: 10.1093/mnras/stx2176.
127
+
128
+
72
129
##### Isotropic (hydro and MHD) and anisotropic thermal conduction (only MHD)
73
130
In the presence of magnetic fields thermal conduction is becoming anisotropic with the flux along
74
131
the local magnetic field direction typically being much stronger than the flux perpendicular to the magnetic field.
@@ -140,6 +197,29 @@ Default value corresponds to the typical value used in literature and goes back
140
197
[^BM82]:
141
198
S. A. Balbus and C. F. McKee, “The evaporation of spherical clouds in a hot gas. III - Suprathermal evaporation,” , vol. 252, pp. 529–552, Jan. 1982, doi: https://doi.org/10.1086/159581
142
199
200
+
#### Viscosity/Momentum diffusion
201
+
202
+
Only isotropic viscosity with a (spatially and temporally) fixed coefficient in code units
203
+
(`code_length`^2/`code_time`) is currently implemented.
204
+
To enable set (in the `<diffusion>` block)
205
+
```
206
+
viscosity = isotropic
207
+
viscosity_coeff = fixed
208
+
mom_diff_coeff_code = 0.25 # fixed coefficent of the kinetmatic viscosity in code units (code_length^2/code_time)
209
+
```
210
+
211
+
#### Resistivity/Ohmic diffusion
212
+
213
+
Only resistivity with a (spatially and temporally) fixed coefficient in code units
214
+
(`code_length`^2/`code_time`)is currently implemented.
215
+
To enable set (in the `<diffusion>` block)
216
+
```
217
+
resistivity = ohmic
218
+
resistivity_coeff = fixed
219
+
ohm_diff_coeff_code = 0.25 # fixed coefficent of the magnetic (ohmic) diffusivity code units (code_length^2/code_time)
220
+
```
221
+
222
+
143
223
### Additional MHD options in `<hydro>` block
144
224
145
225
Parameter: `glmmhd_source` (string)
@@ -243,3 +323,108 @@ d_log_temp_tol = 1e-8 # Tolerance in cooling table between subseque
243
323
Finally, a more comprehensive descriptions of various conventions used for cooling
244
324
functions (and the chosen implementation in AthenaPK) can be found in [Cooling Notes](cooling_notes.md)
245
325
and a notebook comparing various cooling tables (and their conversion) in [cooling/cooling.ipynb](cooling/cooling.ipynb).
326
+
327
+
### Particles
328
+
329
+
#### Tracers
330
+
331
+
Tracer particles can be enabled by setting `tracers/enabled=true` in the input file:
332
+
333
+
```
334
+
<tracers>
335
+
enabled = true
336
+
337
+
initial_seed_method = random_per_block # alternative: user
338
+
initial_num_tracers_per_cell = 0.125
339
+
### Optional arguments
340
+
#initial_rng_seed = INTEGER
341
+
```
342
+
343
+
Two seeding methods are currently supported:
344
+
345
+
-`initial_seed_method=random_per_block`
346
+
- seeds particles at random positions in each block
347
+
- NOTE: the random number generator seed uses the unique block id. Therefore, simulations with the mesh decomposition (mesh and meshblock sizes) are identical independent of the number of MPI ranks used, but if the meshblock size is changed for given mesh (and thus the total number of blocks) the initial state will be different.
348
+
-`initial_num_tracers_per_cell` determines the number of seeded particles per cell
349
+
-`initial_rng_seed` (optional) is used as seed in addition to the block id.
350
+
-`initial_seed_method=user`
351
+
- Calls a problem specific callback function (`ProblemSeedInitialTracers`), see [tracer callback documenation](https://github.com/parthenon-hpc-lab/athenapk/blob/main/docs/pgen.md#tracers).
352
+
353
+
By default, swarm fields are written only to restart files.
354
+
If they are required for "standard" output files (like single precision `hdf5`),
355
+
they need to be added manually to the output block, e.g., (bottom two lines)
356
+
```
357
+
<parthenon/output2>
358
+
file_type = hdf5 # Binary data dump
359
+
variables = prim # variables to be output
360
+
dt = 0.1 # time increment between outputs
361
+
id = prim
362
+
single_precision_output = true
363
+
364
+
swarms = tracers
365
+
tracers_variables = id, x, y, z, rho
366
+
#write_swarm_xdmf=true # uncomment to create an xdmf output (e.g., for Paraview or Visit)
367
+
```
368
+
369
+
Tracers can be read/processed by Paraview or Visit (via the xdmf file)
370
+
or by the `phdf` package shipped with the Parthenon submodule.
371
+
A sample plotting script for the latter might look like
print("Couldn't find module to read Parthenon hdf5 files.")
386
+
387
+
data = phdf.phdf("../run-turbulence/parthenon.prim.00010.phdf")
388
+
tracers = data.GetSwarm("tracers")
389
+
xs = tracers.x
390
+
ys = tracers.y
391
+
zs = tracers.z
392
+
ids = tracers.Get("id")
393
+
394
+
fig = plt.figure()
395
+
ax = fig.add_subplot(projection='3d')
396
+
397
+
ax.scatter(xs, ys, zs,s=1,c=ids)
398
+
```
399
+
resulting in the following image:
400
+
401
+
402
+
403
+
Following "restrictions" apply to the current tracer implementation:
404
+
- Only 3D simulations.
405
+
- Only one advection method (RK2/Heun's method).
406
+
- All primitive fields (`rho`, `pressure`, `vel_x`, `vel_y`, `vel_z`, `B_x`, `B_y`, and `B_x`) are traced by default (independent of whether they're needed or not) in addition to the position (`x`, `y`, `z`) and id (`id`) fields.
407
+
- Default tracer values (such as the primitive fields) are only updated right before writing an output file.
408
+
- Ids are only unique if tracers are seeded at the beginning at the simulations and no new tracer particles are added dynamically while the simulation is running.
409
+
410
+
Please get in touch, if you interested in running simulations that require lifting one (or more) of those restrictions.
411
+
412
+
## Boundary conditions
413
+
414
+
In addition to enrolling custom boundary conditions, three general options are currently supported by default:
415
+
416
+
-`periodic`: Work without restrictions (implemented in upstream Parthenon)
417
+
-`outflow`: Last layer of active cells is copied into ghost cells. Work (technicaly) without restrictions (implemented in upstream Parthenon). Physically care must be taken when using outflow conditions in subsonic flows (as characteristics can travel into the domain, which is not properly handled by copying the active cell data).
418
+
It is recommended to implement custom, problem specific boundary conditions for subsonic outflows, see large body of engineering literature.
419
+
For supersonic outflows there is no issue as the flow itself is faster than the fastest characteristic.
420
+
-`reflecting`: Work only for cell-centered *hydro* fluid variables
421
+
(implemented in AthenaPK currently without support for particles or `Metadata::Fine` fields,
422
+
where the latter are currently not being used in AthenaPK).
423
+
Note, that we specifically do *not* provide "reflecting" boundary conditions for the magnetic
424
+
fields as these kind of configurations are typically highly problem dependent and, thus,
425
+
require manual treatment.
426
+
427
+
They are being used by specifying the name as parameter value for the
428
+
`ix1_bc`, `ox1_bc`, `ix2_bc`, `ox2_bc`, `ix3_bc`, and `ox3_bc` paraemters
429
+
in the `<parthenon/mesh>` block of the input file for the inner and outer
430
+
(say left and right in x1 direction) in each coordinate direction, respectively.
0 commit comments