Skip to content

Add code for modeling a radiation pressure dominated disk#14

Open
soumide1102 wants to merge 14 commits into
lanl:masterfrom
soumide1102:rad_press_torus1
Open

Add code for modeling a radiation pressure dominated disk#14
soumide1102 wants to merge 14 commits into
lanl:masterfrom
soumide1102:rad_press_torus1

Conversation

@soumide1102

Copy link
Copy Markdown
Collaborator

This PR adds code for setting up a radiation pressure dominated disk.

  • For a gamma disk, the code now can either set up a radiation pressure dominated or a gas pressure dominated case. For doing that the user can input -gammagaspress or gammaradpress. -gammagaspress will follow the old code path for gamma. -gammaradpress will follow the radiation pressure code path added in this PR.
  • For the radiation-pressure dominated case, rho and u in problem.c is calculated using an inputted entropy.
  • Also modified the other files in nublight to incorporate the radiation pressure dominated condition.
  • The kappa parameter saved in build.py is now corrected to kappa_eos

In calculating rho for the radiation pressure dominated case, there is a multiplication factor a that is currently not included in problem.c . Should hdf5_to_dict be the place to factor that in when loading the results and converting to the desired units?

@soumide1102 soumide1102 requested a review from jonahm-LANL March 1, 2021 20:56
Comment thread core/fixup.c Outdated
Comment on lines +232 to +234
#if ELECTRONS && EOS == EOS_TYPE_GAMMA_GASPRESS
// Reset entropy after floors
pv[KTOT] = EOS_Gamma_entropy_rho0_u(pv[RHO], pv[UU]);
pv[KTOT] = EOS_Gamma_Gaspress_entropy_rho0_u(pv[RHO], pv[UU]);

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I also wanted to check with you what this piece of code is doing. In the current modifications, I changed the name of the functions to the updated names specifying gas pressure. If this is going into evolving the simulation, should there also be another if statement for the radiation pressure dominated case as well?

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This does go into evolving the simulation, yes. But this code path will never be reached. The ELECTRONS machinery is backwards compatibility with ebhlight's two-temperature gas for AGN disks. I would add another if case, but it's never going to come up.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

With EOS_TYPE_GAMMA having two sub-cases, this can go to EOS_entropy_rho0_u or something

@Yurlungur Yurlungur left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Overall I like the strategy here... But something seems off as far as the code paths. You only replaced the entropy equation for the ideal gas equation of state, but you're treating it as an entirely different equation of state everywhere in the code. I think this is only working because you set all 3 macros to resolve to the same integer.

I would suggest a new structure. Something like:

  1. The ideal gas eqation of state is triggered for EOS == EOS_TYPE_GAMMA
  2. The gamma EOS has two new options GAMMA_GASPRESS and GAMMA_RADPRESS or something.
  3. Those new macros are set as a new cparm in the build fie. We probably also need to set a new default value in script/configure.py

Comment thread core/decs.h Outdated
Comment on lines +199 to +200
#define EOS_TYPE_GAMMA_GASPRESS (0)
#define EOS_TYPE_GAMMA_RADPRESS (0)

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These all can't be zero. That means EOS_TYPE_GAMMA_GASPRESS and EOS_TYPEGAMMA_RADPRESS will resolve to the same number. That's probably not your intent. Each EOS type needs to resolve to a different number in the macro.

Comment thread core/main.c Outdated
fprintf(stdout, " * IDEAL GAS EOS *\n");
#endif
#if EOS == EOS_TYPE_GAMMA_GASPRESS
fprintf(stdout, " * GAS-PRESSURE DOMINATED EOS *\n");

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

GAS-PRESSURE DOMINATED IDEAL GAS EOS

Comment thread core/main.c Outdated
fprintf(stdout, " * GAS-PRESSURE DOMINATED EOS *\n");
#endif
#if EOS == EOS_TYPE_GAMMA_RADPRESS
fprintf(stdout, " * RADIATION-PRESSURE DOMINATED EOS *\n");

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

RADIATION-PRESSURE DOMINATED IDEAL GAS EOS

Comment thread prob/torus_cbc/build.py Outdated
Rout_vis = Rout
GAMMA = 13./9.
KAPPA = 1.e-3
GAMMA = 4./3.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should be if-cased on whether you're radiation or gas pressure dominated. 4./3. is appropriate for radiation pressure, not gas pressure.

Comment thread prob/torus_cbc/build.py Outdated
elif DO_GAMMA_RADPRESS:
EOS_TYPE = "EOS_TYPE_GAMMA_RADPRESS"
elif DO_GAMMA_GASPRESS:
EOS_TYPE = "EOS_TYPE_GAMMA_GASPRESS"

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

now there should be a throw if it's not one a valid EOS type. Let's catch that in case someone somehow doesn't set one.

@soumide1102

Copy link
Copy Markdown
Collaborator Author

@Yurlungur Does this look like what you were suggesting? I will then update the other scripts accordingly.

Comment thread prob/torus_cbc/build.py

# EOS
bhl.config.set_cparm("EOS_TYPE", EOS_TYPE)
bhl.config.set_cparm("EOS", EOS_TYPE)

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So one small change I was thinking of: Here the name of the parameter is EOS which saves which type of EoS is going to be used. Also, around line 100, there is

if RELTABLE:
    bhl.report_var('EOS','RELTABLE')
elif GAMTABLE:
    bhl.report_var('EOS','GAMTABLE')
else:
    bhl.report_var('EOS','GAMMA')

where the parameter name is EOS. I was wondering whether the latter parameter can be changed to something like EOS_VAL to keep the two parameter names different?

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

no those should be the same. report_var doesn't do anything internal in the code. It prints a choice of variable to the build output.

Comment thread core/decs.h
Comment on lines +204 to +206
#define GASPRESS (1)
#define RADPRESS (2)
#elif EOS == EOS_TYPE_POLYTROPE

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes. This is what I was thinking. :)

Comment thread core/eos.c Outdated
Comment on lines +162 to +166
#if EOS == EOS_TYPE_GAMMA_GASPRESS
ent = EOS_Gamma_Gaspress_entropy_rho0_u(rho, u);
#elif EOS == EOS_TYPE_GAMMA_RADPRESS
ent = EOS_Gamma_Radpress_entropy_rho0_u(rho, u);

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should now be sub-cases of EOS_TYPEGAMMA

Comment thread core/fixup.c Outdated
Comment on lines +232 to +234
#if ELECTRONS && EOS == EOS_TYPE_GAMMA_GASPRESS
// Reset entropy after floors
pv[KTOT] = EOS_Gamma_entropy_rho0_u(pv[RHO], pv[UU]);
pv[KTOT] = EOS_Gamma_Gaspress_entropy_rho0_u(pv[RHO], pv[UU]);

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

With EOS_TYPE_GAMMA having two sub-cases, this can go to EOS_entropy_rho0_u or something

Comment thread core/io.c Outdated
#endif // METRIC

#if EOS == EOS_TYPE_GAMMA || GAMMA_FALLBACK
#if EOS == EOS_TYPE_GAMMA || EOS_TYPE_GAMMA_GASPRESS || EOS_TYPE_GAMMA_RADPRESS || GAMMA_FALLBACK

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can go back to just EOS == EOS_TYPE_GAMMA

Comment thread core/io.c Outdated

// EOS
#if EOS == EOS_TYPE_GAMMA || GAMMA_FALLBACK
#if EOS == EOS_TYPE_GAMMA || EOS_TYPE_GAMMA_GASPRESS || EOS_TYPE_GAMMA_RADPRESS || GAMMA_FALLBACK

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can go back to just EOS == EOS_TYPE_GAMMA

Comment thread core/io.c Outdated
READ_HDR(cour, TYPE_DBL);

#if EOS == EOS_TYPE_GAMMA || GAMMA_FALLBACK
#if EOS == EOS_TYPE_GAMMA || EOS_TYPE_GAMMA_GASPRESS || EOS_TYPE_GAMMA_RADPRESS || GAMMA_FALLBACK

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can go back to just EOS == EOS_TYPE_GAMMA

@Yurlungur

Copy link
Copy Markdown
Collaborator

@soumide1102 this is looking like it's moving along... although I see tests are failing. How close are you to being ready for another review?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants