Question about code units and physical units (hydrodynamics simulations)

Question

I'm working on a code that implements smoothed particle hydrodynamics (SPH) method for solving the equations of magnetohydrodynamics (MHD) with self-gravity.

In research papers regarding existing software the authors sometimes mention storing computational values in "code units" as opposed to physical units. Also, some test problems are given with initial conditions without units like "M = 1" (total mass; I'm assuming here that this means 1 [unit-of-mass] in code units).

Specifically, in gravitational hydrodynamics the gravitational constant may be set simply as 1; in MHD the permeability of free space constant can be set as 1.

From my understanding this means that all physical values need to be scaled accordingly to give us the same initial conditions represented in code units.

For example, if we have a sphere of gas of radius R_phys = 1e16 [cm] we would have to scale it to get

R_code = a * R_phys = x [length],

where a is some constant relating physical units of length and code units.

But this doesn't require any changes in the physical equations, right? We're supposed to simply define the gravitational constant and permeability of free space as unity and all we have to worry about is scaling of all physical values to correctly represent a specific test problem given in physical units? (Oh, and of course, we would need to scale code units back to physical if we want to compare the simulation results with some observational data).

Let's say that I don't care about taking a specific test problem (given in physical units) and I just set the gravitational constant and permeability of free space as unity. Let's say that I set domain range as 1 [length] and total mass as 1 [mass] (which would give certain value of gravitational energy for the gas cloud). And from there I set initial temperature/pressure, strength of magnetic field and angular velocity just by estimating the ratios of corresponding energies to the gravitational energy (to relate it to physical problems with known ranges of those "energy" parameters).

Note: from my understanding, storing physical values in code units is not important since modern compilers can handle operating with values of great range (specifically for my simulations, particle density can be of order $10^{-17}$ g/cm^3 while, say, magnetic field is ~ $10^{-6}$ gauss, etc.); however I'm not sure about this, and overall it would be better to have all values stored in code units.

Answer 1

But this doesn't require any changes in the physics equations, right?

Correct, but you should convince yourself this is true by choosing the appropriate scaling (I like to use time, length & density as my base 3, as I mention here in a related question) and rearranging the variables in the system of equations to (a) determine the scaling factors for the other variables/units and (b) show that the equations don't change form (just variable names/labels).

You do want to be careful about using $G=1$, because I believe that does require specific choices of units that may be inconsistent with other choices you might want to make. To be sure, you want to set the scale factors for density, velocity, etc. and then reduce $G$ to its dimensionless form by canceling out the units (this is also true for all dimensionful constants).

from my understanding, storing physical values in code units is not important since modern compilers can handle operating with values of great range...

Sorta. Sure you can do 1e-17 * 1e-6 in code, but you're likely going to lose precision as compared to say multiplying 1e-2 * 1e-3 (assuming appropriate scaling factors here to make the two products physically equivalent). This is more prevalent when you're computing the momentum of a gas cloud with $\rho\sim10^{-17}$ g/cc and $v\sim10^4$ km/s, rather than products near unity, but the impact can be seen elsewhere. If you haven't already, you do want to read David Goldberg's What Every Computer Scientist Should Know About Floating-Point Arithmetic.