Cosmo dynamical de

[astrodiegog]

Include ability to calculate the expansion history of cosmological simulations with a time-evolving dark energy equation of state

• Create new parameter expected as wDE_file which points to a file with z vs wDE(z) entries
• Create macro DYNAMICAL_DE_TABLE so that compiler knows whether to allocate certain variables
• Create routine that integrates wDE(z) to find rhoDE(z)/rhoDE(z=0) using midpoint rule
• Linearly interpolate between rhoDE(z)/rhoDE(z=0) entries for arbitrary redshifts
• In verifying the expansion history file output, realized the the expansion history was updated before output : Cosmo.dt += dt_secs occurred before Cosmo.Write_Expansion_History_Entry(). Updated so that expansion history entry occurs before updating the time
• Compilation with DYNAMICAL_DE_TABLE overrides w0, wa parameter file entries

[astrodiegog]

pre-commit.ci autofix

[astrodiegog]

@brantr Comment whenever you have the time!

[evaneschneider]

I'll take a more thorough look at this tomorrow (we have a hack day in the afternoon), but a couple quick comments:
-Since the first CAAR project, we've been trying hard to decrease macro use, because they make it much harder to test code, and are easy to forget to define, or define incorrectly. That doesn't mean we won't allow new compiler macros, but I am trying to encourage input parameter usage as opposed to macros whenever possible. So any new macros should be justified (explain in the PR why you need a macro to accomplish this). That way, if we can help you rewrite the code to NOT use a macro, we will have the opportunity to do so.
-Now that we have read the docs integrated into the codebase, I am enforcing documentation with PRs that include new code capabilities. So this PR should be modified to include, at a minimum, a page (probably in the physics subdirectory) that explains the dynamical dark energy capability AND how it is implemented (the new macro needs to be documented somewhere, there should be a minimum description of the contents of the expected DE file...). I know this may seem silly since as of right now the cosmology portion of Cholla is virtually undocumented, but we have to start somewhere.

[mabruzzo]

(Feel free to ignore me) I took a quick look, and at a glance, it looks like it would be relatively straightforward to update this PR in order to get rid of the DYNAMICAL_DE_TABLE table.

The simplest way to get rid of the parameter is to update the code so that the Parse_Params function always checks (when the COSMOLOGY macro is defined) to see whether the wDE_file parameter is specified. If it's specified, then you will use. If it isn't specified, you could just store an empty string within Parameters::wDE_file.

I included my recommendations for how you might do this down below. (I'm pretty confident it will "just work" if you go through and simply accept all the suggestions, but I haven't actually tested the code).

I do have a couple stylistic suggestions - but I'll hold off on those until after you decide on whether or not to remove the macro.

EDIT: OOPs! I realized that I forgot mention that this looks fantastic! It looks like a wonderful addition to Cholla!

[evaneschneider]

On a related note, this PR fails the cosmology system test, because it changes the default cosmology build to one with dynamical dark energy. That's good! (Presumably, adding dynamical dark energy should change the results of a cosmological simulation!) But, you probably don't actually want to change the default cosmology system test. Simply removing the new macro from the cosmology build should fix this, but that's one of the downsides to approaching this change with a macro - if you wanted to add a new system test for dynamical dark energy (not a bad idea!), we'd have to add a whole new build type to the testing matrix. @mabruzzo made some suggestions for how you could accomplish this with an input parameter rather than a macro. Let us know if anything is confusing!

[mabruzzo]

Now that we have read the docs integrated into the codebase, I am enforcing documentation with PRs that include new code capabilities. So this PR should be modified to include, at a minimum, a page (probably in the physics subdirectory) that explains the dynamical dark energy capability AND how it is implemented (the new macro needs to be documented somewhere, there should be a minimum description of the contents of the expected DE file...). I know this may seem silly since as of right now the cosmology portion of Cholla is virtually undocumented, but we have to start somewhere.

Just in case it helps, this page (here) should help you get started with adding docs.

Please reach out if have any questions or want help (with writing docs or the suggested changes); I'm more than happy to help!

[astrodiegog]

pre-commit.ci autofix

[astrodiegog]

The Cosmology would then track a pointer to this class and call the methods of that class.

I am unsure how best to allocate the memory for this pointer, so I've just saved the method as an object instead of a pointer to this object.

I believe I addressed other comments related to adding some documentation, creating a separate class for the dynamical dark energy implementation, and using std::vector<float> instead of just pointers to arrays.

I am happy to take further comments when people have time!

[mabruzzo]

Hi Diego,

This looks exceptionally! I have left a few comments. You can totally ignore all the ones that aren't marked IMPORTANT (the IMPORTANT comments only ask for very minor changes)

EDIT: actually, I meant to ask why you have dynamical_dark_energy.h and tabulated_dynamicalDE_EoS.cpp (rather than dynamical_dark_energy.h and dynamical_dark_energy.cpp OR tabulated_dynamicalDE_EoS.h and tabulated_dynamicalDE_EoS.cpp). Was that intentional?

[mabruzzo]

(THIS TOTALLY ISN'T NECESSARY), but since you've gone to the trouble of shifting this logic into a nice class, it might be nice to add a constructor:

Suggested change

	/*! Load redshift and dark energy equation of state table z, wDE(z), only called once to setup dynamical DE case */
	void Setup_DynamicalDE_EquationOfState_(struct Parameters *P);
	/*! Calculate dark energy density normalized to z=0, populate dynamicalDE_table_density */
	void Set_DynamicalDE_Density();
	private:
	/*! Load redshift and dark energy equation of state table z, wDE(z), only called once to setup dynamical DE case */
	void Setup_DynamicalDE_EquationOfState_(struct Parameters *P);
	/*! Calculate dark energy density normalized to z=0, populate dynamicalDE_table_density */
	void Set_DynamicalDE_Density();
	public:
	explicit TabulatedDynamicalDarkEnergyEoS(Parameters *P)
	{
	Setup_DynamicalDE_EquationOfState_(P);
	Set_DynamicalDE_Density();
	}

In case you didn't know:

• we can pass Parameters *P rather than struct Parameters *P (the latter is only necessary in pure C)
• it's generally good practice to annotate a constructor that accepts a single argument with the explicit keyword. If you don't do this C++ is allowed to perform implicit casts (while this sounds convenient, it is actually a bad thing because it can lead to very unexpected behavior)

[mabruzzo]

If you decide to make this change, then you would need to replace

tab_dynamicalDE_EoS.Setup_DynamicalDE_EquationOfState_(P);
tab_dynamicalDE_EoS.Set_DynamicalDE_Density();

within cosmology.cpp with tab_dynamicalDE_EoS = TabulatedDynamicalDarkEnergyEoS(P);

[astrodiegog]

Having a constructor makes the code look much cleaner, but I am having trouble integrating this. When I compile, I receive the following errors:

src/cosmology/cosmology.cpp: In constructor ‘Cosmology::Cosmology()’:
src/cosmology/cosmology.cpp:7:26: error: no matching function for call to ‘TabulatedDynamicalDarkEnergyEoS::TabulatedDynamicalDarkEnergyEoS()’
    7 | Cosmology::Cosmology(void) {}
      |                          ^
In file included from src/cosmology/../cosmology/cosmology.h:13,
                 from src/cosmology/cosmology.cpp:3:
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:22:12: note: candidate: ‘TabulatedDynamicalDarkEnergyEoS::TabulatedDynamicalDarkEnergyEoS(Parameters*)’
   22 |   explicit TabulatedDynamicalDarkEnergyEoS(struct Parameters *P);
      |            ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:22:12: note:   candidate expects 1 argument, 0 provided
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:10:7: note: candidate: ‘TabulatedDynamicalDarkEnergyEoS::TabulatedDynamicalDarkEnergyEoS(const TabulatedDynamicalDarkEnergyEoS&)’
   10 | class TabulatedDynamicalDarkEnergyEoS
      |       ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:10:7: note:   candidate expects 1 argument, 0 provided
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:10:7: note: candidate: ‘TabulatedDynamicalDarkEnergyEoS::TabulatedDynamicalDarkEnergyEoS(TabulatedDynamicalDarkEnergyEoS&&)’
src/cosmology/../cosmology/tabulated_dynamicalDE_EoS.h:10:7: note:   candidate expects 1 argument, 0 provided

I think this has to do with the definition TabulatedDynamicalDarkEnergyEoS tab_dynamicalDE_EoS; in the Cosmology class being created without using its constructor. From what I found online, people seem to recommend storing a pointer instead of the object itself.

This would also tie in with your most recent comment about using std::unique_ptr to store a pointer to this table !

[mabruzzo]

You are exactly right! This pertains to default constructors.

I can provide you with a little more context:

• the error is pointing to the definition of Cosmology's default constructor. Since Cosmology's default constructor doesn't invoke any particular constructor for its tab_dynamicalDE_EoS data member,¹ the compiler will try to invoke that data member's default constructor.
• in the code's current state (without the change I suggested here), the fact that you haven't defined any constructor for TabulatedDynamicalDarkEnergyEoS means that the compiler will implicitly define a default constructor that simply calls the default constructor of each data member in TabulatedDynamicalDarkEnergyEoS
• when you implement the change I suggested here, the fact that you are implementing any constructor for TabulatedDynamicalDarkEnergyEoS means that the compiler stops implementing a default constructor TabulatedDynamicalDarkEnergyEoS

There are several ways to resolve this issues

1. Simple approach: just add a default-constructor to TabulatedDynamicalDarkEnergyEoS. This might look like:

    public:
     TabulatedDynamicalDarkEnergyEoS() = default;
   
     explicit TabulatedDynamicalDarkEnergyEoS(Parameters *P)
     {
       Setup_DynamicalDE_EquationOfState_(P);
       Set_DynamicalDE_Density();
     }

   There are multiple other ways to do this.²

2. As you already described, you can convert Cosmology::tab_dynamicalDE_EoS so that it's a std::unique_ptr<TabulatedDynamicalDarkEnergyEoS>.

   • This is what I would recommend.
   • If you have trouble using std::unique_ptr, it may be easier to first write a commit that uses regular pointers (i.e. a TabulatedDynamicalDarkEnergyEoS*) and then write a subsequent commit with std::unique_ptr. (Also feel free to ask me for help).
   • The reason this works is that std::unique_ptr has a default constructor -- by default it stores a null pointer.³

3. You can convert Cosmology::tab_dynamicalDE_EoS so that it's a std::optional<TabulatedDynamicalDarkEnergyEoS>.

   • this theoretically provides certain performance advantages over using a std::unique_ptr, but they aren't significant in the current context.⁴
   • since you are still familiarizing yourself with C++, it probably makes more sense for you to use a std::unique_ptr

4. You can convert Cosmology::tab_dynamicalDE_EoS so that it's a container (e.g. a std::vector) of TabulatedDynamicalDarkEnergyEoS instances. This doesn't make any real sense in the current context.

Personally, I generally view solutions like options 2 and 3 as more elegant than option 1:

• If you have a default constructor (option 1), its possible to have an instance of TabulatedDynamicalDarkEnergyEoS that is empty. In this case, its the programmer's has to know whether the table is empty or not to decide whether they can call the Get_DynamicalDE_Density_from_a method.
• If you don't have a default constructor (options 2 or 3) then it's impossible to have an empty table
• practically speaker, the distinction isn't very significant since the class is so simple

Footnotes

1. ASIDE when defining a class's constructor, you can specify which of the constructors to call for creating a data member using member initializer lists. ↩

2. The = default; syntax that I used is a way to explicitly tell the compiler to define the same default constructor that it had been implicitly defining before this change. ↩

3. If you use a regular pointer, it doesn't have any constructor. It will be initialized in an indeterminate state (just like in C). ↩

4. For some context, a std::unique_ptr or std::optional can both be used to signal that it wraps a value that may or may not exist. If the value does exist, std::unique_ptr holds a pointer to the memory containing that value whereas a std::optional directly holds the memory containing the value. Thus, using a std::optional can be faster since there are fewer memory dereferences. In the way that Cosmology::tab_dynamicalDE_EoS is currently used (and the overhead of branch-prediction), I doubt this performance differential really matters. ↩

[astrodiegog]

Gotcha ! This makes sense to me. Thank you for taking the time to provide some detailed context, and the different ways to go about solving this issue.

I opted for option 2 (using std::unique_ptr)! I just pushed these changes

[mabruzzo]

The Cosmology would then track a pointer to this class and call the methods of that class.

I am unsure how best to allocate the memory for this pointer, so I've just saved the method as an object instead of a pointer to this object.

I believe I addressed other comments related to adding some documentation, creating a separate class for the dynamical dark energy implementation, and using std::vector<float> instead of just pointers to arrays.

I am happy to take further comments when people have time!

Just in case your curious, the way I would handle the pointer allocation stuff is I would make a unique pointer. This would involve:

• including <memory> at the top of comology.h
• replacing TabulatedDynamicalDarkEnergyEoS tab_dynamicalDE_EoS; in cosmology.h with std::unique_ptr<TabulatedDynamicalDarkEnergyEoS> tab_dynamicalDE_EoS;
• replacing

    tab_dynamicalDE_EoS.Setup_DynamicalDE_EquationOfState_(P);
    tab_dynamicalDE_EoS.Set_DynamicalDE_Density();

  with tab_dynamicalDE_EoS = std::unique_ptr<TabulatedDynamicalDarkEnergyEoS>( new TabulatedDynamicalDarkEnergyEoS(P));
  (this assumes that you added a constructor like I mentioned int the comments)

You would be able to check if the underlying pointer is empty by calling tab_dynamicalDE_EoS.get() == nullptr. And you could call methods on this pointer by writing tab_dynamicalDE_EoS->Get_DynamicalDE_Density_from_a(a);

(To be clear, you definitely don't have to do this -- I just wanted to explain how it would be accomplished)

[astrodiegog]

Was that intentional?

I was flipping between both and forgot to decide ! I have changed dynamical_dark_energy.h to tabulated_dynamicalDE_EoS.h

[astrodiegog]

pre-commit.ci autofix

[astrodiegog]

Would like to work on storing TabulatedDynamicalDarkEnergyEoS as pointer in Cosmology in next couple days!

[mabruzzo]

@astrodiegog, Now that you have converted TabulatedDynamicalDarkEnergyEoS into its own class, it would be great if we could write a simple unit test.

You've already gone above and beyond in this PR, and learning how to use GoogleTest can be a little involved. So if you give me a sample file that this class reads, I would be happy to do that myself (I can open a new PR that would be reviewed after this one is merged)

[astrodiegog]

pre-commit.ci autofix

[astrodiegog]

You've already gone above and beyond in this PR, and learning how to use GoogleTest can be a little involved. So if you give me a sample file that this class reads, I would be happy to do that myself (I can open a new PR that would be reviewed after this one is merged)

Thanks ! I just sent a sample file via Slack