enzo-project / enzo-e

A version of Enzo designed for exascale and built on charm++.
Other
29 stars 35 forks source link

Inconsistent cosmology units before/after PR #234 #243

Closed WillHicks96 closed 2 years ago

WillHicks96 commented 2 years ago

I noticed this when I was preparing PR #242.

If I rollback to the version of main just before PR #234, run test_cosmo-dd.in, and print out all of the cosmology units during cycle 0, I get these values:

cosmology->density_units()  = 2.76e-24
cosmology->length_units()   = 1.32e+23
cosmology->time_units()     = 6.57e+14
cosmology->velocity_units() = 2.01e+08

If I do the same for the current version of main, I get the following:

cosmology->density_units()  = 3.96e-17
cosmology->length_units()   = 1.34e+23
cosmology->time_units()     = 6.57e+14
cosmology->velocity_units() = 2.01e+09 

I went through the math outlined in EnzoPhysicsCosmology::density_units(), and it looks like an extra factor of G was picked up in the denominator at some point (I think it should be G in the denominator, not G^2), which can account for the factor of ~1e7ish difference in density_units(). Fixing this in the current version of main gives the correct value of density_units()=2.76e-24.

velocity_units() is still off by a factor of 10 though, and I don't quite follow where the final equation returned by EnzoPhysicsCosmology::velocity_units() comes from. Is there any reason why we can't just set velocity_units() = length_units()/time_units()? The numbers seem to work out in this case.

mabruzzo commented 2 years ago

I went through the math outlined in EnzoPhysicsCosmology::density_units(), and it looks like an extra factor of G was picked up in the denominator at some point (I think it should be G in the denominator, not G^2), which can account for the factor of ~1e7ish difference in density_units(). Fixing this in the current version of main gives the correct value of density_units()=2.76e-24.

Sounds sensible to me. Somebody more knowledgeable than me should comment on this.

velocity_units() is still off by a factor of 10 though, and I don't quite follow where the final equation returned by EnzoPhysicsCosmology::velocity_units() comes from. Is there any reason why we can't just set velocity_units() = length_units()/time_units()? The numbers seem to work out in this case.

I think you're absolutely right that we should set velocity_units = length_units()/time_units()

stefanarridge commented 2 years ago

cosmology->density_units() = 2.76e-24 cosmology->length_units() = 1.32e+23 cosmology->time_units() = 6.57e+14 cosmology->velocity_units() = 2.01e+08


If I do the same for the current version of `main`, I get the following:

cosmology->density_units() = 3.96e-17 cosmology->length_units() = 1.34e+23 cosmology->time_units() = 6.57e+14 cosmology->velocity_units() = 2.01e+09



I went through the math outlined in `EnzoPhysicsCosmology::density_units()`, and it looks like an extra factor of G was picked up in the denominator at some point (I think it should be G in the denominator, not G^2), which can account for the factor of ~1e7ish difference in `density_units()`. Fixing this in the current version of `main` gives the correct value of `density_units()=2.76e-24`.

Yes, you're right. The extra factor of G was a mistake on my part.

velocity_units() is still off by a factor of 10 though, and I don't quite follow where the final equation returned by EnzoPhysicsCosmology::velocity_units() comes from. Is there any reason why we can't just set velocity_units() = length_units()/time_units()? The numbers seem to work out in this case.

The velocity units are (1 + z_i) / (1+z) * length_unit() / time_unit(). z_i is the initial redshift, and z is the current redshift. (Edit: I just realised I made a mistake, there should be a sqrt(1+initial_redshift) instead of a (1+initial_redshift)) I put a comment in src/Enzo/enzo_EnzoPhysicsCosmology.hpp explaining the reasoning as I understood it. @gregbryan is probably a better person to ask regarding why this choice was made in Enzo.

Also, I find it slightly strange that the length unit is slightly different since I didn't make any changes to this.

WillHicks96 commented 2 years ago

I just realised I made a mistake, there should be a sqrt(1+initial_redshift) instead of a (1+initial_redshift)

@stefanarridge Ahh okay, that explains the factor of 10 difference in velocity_units -- The initial redshift is 99 for test_cosmo