From b7394c09209ab25df6cb05570ae04cc3a11390b8 Mon Sep 17 00:00:00 2001 From: Ben Wibking Date: Mon, 4 Nov 2024 17:01:47 -0500 Subject: [PATCH 1/2] add sqrt(4pi) to cgs B-field units --- src/units.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/units.hpp b/src/units.hpp index 3290e87e..7579354d 100644 --- a/src/units.hpp +++ b/src/units.hpp @@ -90,7 +90,7 @@ class Units { return code_energy_cgs() / kev_cgs * code_length_cgs() * code_length_cgs(); } parthenon::Real code_magnetic_cgs() const { - return sqrt(code_mass_cgs()) / sqrt(code_length_cgs()) / code_time_cgs(); + return std::sqrt(4.0 * M_PI) * sqrt(code_mass_cgs()) / sqrt(code_length_cgs()) / code_time_cgs(); } // Physical Constants in code units From 34be9049dcfa2b5e9008ee82079aed6952b6daf2 Mon Sep 17 00:00:00 2001 From: Philipp Grete Date: Fri, 8 Nov 2024 16:52:49 +0100 Subject: [PATCH 2/2] Add unit doc --- CHANGELOG.md | 1 + README.md | 17 ++------------- docs/README.md | 1 + docs/input.md | 4 ++++ docs/units.md | 58 ++++++++++++++++++++++++++++++++++++++++++++++++++ src/units.hpp | 3 ++- 6 files changed, 68 insertions(+), 16 deletions(-) create mode 100644 docs/units.md diff --git a/CHANGELOG.md b/CHANGELOG.md index 152bd814..c95a86ae 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,7 @@ - [[PR 1]](https://github.com/parthenon-hpc-lab/athenapk/pull/1) Add isotropic thermal conduction and RKL2 supertimestepping ### Changed (changing behavior/API/variables/...) +- [[PR 122]](https://github.com/parthenon-hpc-lab/athenapk/pull/122) Fixed sqrt(4pi) factor in CGS Gauss unit and add unit doc - [[PR 119]](https://github.com/parthenon-hpc-lab/athenapk/pull/119) Fixed Athena++ paper test case for KHI pgen. Added turbulence pgen doc. - [[PR 97]](https://github.com/parthenon-hpc-lab/athenapk/pull/97) Fixed Schure cooling curve. Removed SD one. Added description of cooling function conventions. - [[PR 84]](https://github.com/parthenon-hpc-lab/athenapk/pull/84) Bump Parthenon to latest develop (2024-02-15) diff --git a/README.md b/README.md index fd66e824..882d9625 100644 --- a/README.md +++ b/README.md @@ -131,21 +131,8 @@ the `file_type = hdf5` format, see [VisIt](https://wci.llnl.gov/simulation/computer-codes/visit/). In ParaView, select the "XDMF Reader" when prompted. -2. With [yt](https://yt-project.org/) -- though currently through a custom frontend -that is not yet part of the main yt branch and, thus, has to be installed manually, e.g., -as follows: -```bash -cd ~/src # or any other folder of choice -git clone https://github.com/forrestglines/yt.git -cd yt -git checkout parthenon-frontend - -# If you're using conda or virtualenv -pip install -e . -# OR alternatively, if you using the plain Python environment -pip install --user -e . -``` -Afterwards, `*.phdf` files can be read as usual with `yt.load()`. +2. With [yt](https://yt-project.org/) +As of versions >=4.4 `*.phdf` files can be read as usual with `yt.load()`. 3. Using [Ascent](https://github.com/Alpine-DAV/ascent) (for in situ visualization and analysis). This requires Ascent to be installed/available at compile time of AthenaPK. diff --git a/docs/README.md b/docs/README.md index 018d4b68..145cc0f8 100644 --- a/docs/README.md +++ b/docs/README.md @@ -13,6 +13,7 @@ The documentation currently includes - [Notebooks to calculate cooling tables from literature](cooling) - [Brief notes on developing code for AthenaPK](development.md) - [How to add a custom/user problem generator](pgen.md) +- [Units](units.md) - Detailed descriptions of more complex problem generators - [Galaxy Cluster and Cluster-like Problem Setup](cluster.md) - [Driven turbulence](turbulence.md) diff --git a/docs/input.md b/docs/input.md index d5e48d4f..73f2b23b 100644 --- a/docs/input.md +++ b/docs/input.md @@ -67,6 +67,10 @@ To control the floors, following parameters can be set in the `` block: *Note* the pressure floor will take precedence over the temperature floor in the conserved to primitive conversion if both are defined. +#### Units + +See(here)[units.md]. + #### Diffusive processes Diffusive processes in AthenaPK can be configured in the `` block of the input file. diff --git a/docs/units.md b/docs/units.md new file mode 100644 index 00000000..b26ccf6e --- /dev/null +++ b/docs/units.md @@ -0,0 +1,58 @@ +# AthenaPK units + +## General unit system + +Internally, all calculations are done in "code units" and there are no conversions between +code and physical units during runtime (with excetion like the temperature when cooling is +being used). +Therefore, in general no units need to be prescribed to run a simulation. + +If units are required (e.g., if cooling is used and, thus, a conversion between internal energy +in code units and physical temperature is required) they are configured in the input block +as follows: + +``` + +code_length_cgs = 3.085677580962325e+24 # 1 Mpc in cm +code_mass_cgs = 1.98841586e+47 # 1e14 Msun in g +code_time_cgs = 3.15576e+16 # 1 Gyr in s +``` + +This information will also be used by postprocessing tools (like yt) to convert between +code units and a physical unit system (like cgs). + +Moreover, internally a set of factors from code to cgs units are available to process conversions +if required (e.g., from the input file). + +For example, for an input parameter (in the input file) like + +``` + +r0_cgs = 3.085677580962325e+20 # 100 pc +``` + +the conversion should happen in the problem generator lik + +```c++ + r_cloud = pin->GetReal("problem/cloud", "r0_cgs") / units.code_length_cgs(); +``` + +so that the resulting quantity is internally in code units (here code length). + +It is highly recommended to be *very* explicit/specific about units everywhere (as it is +a common source of confusion) like adding the `_cgs` suffix to the parameter in the +input file above. + +## Magnetic units + +Internally, AthenaPK (and almost all MHD codes) use +[Heaviside-Lorentz units](https://en.wikipedia.org/wiki/Heaviside%E2%80%93Lorentz_units), +where the magnetic field is transformed from $B \rightarrow B / \sqrt{4 \pi}$. +(See also the note in the +[Castro documentation](https://amrex-astro.github.io/Castro/docs/mhd.html) about this.) + +So when converting from CGS-Gaussian units to code units, it is necessary to divide +by $\sqrt{4 \pi}$ (in addition to the base dimensional factors). +This is automatically handled by the `units.code_magnetic_cgs()` factors. + + diff --git a/src/units.hpp b/src/units.hpp index 7579354d..fc6a2412 100644 --- a/src/units.hpp +++ b/src/units.hpp @@ -90,7 +90,8 @@ class Units { return code_energy_cgs() / kev_cgs * code_length_cgs() * code_length_cgs(); } parthenon::Real code_magnetic_cgs() const { - return std::sqrt(4.0 * M_PI) * sqrt(code_mass_cgs()) / sqrt(code_length_cgs()) / code_time_cgs(); + return std::sqrt(4.0 * M_PI) * sqrt(code_mass_cgs()) / sqrt(code_length_cgs()) / + code_time_cgs(); } // Physical Constants in code units