Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

calculating Mass accretion rate for SN post-processing #211

Merged
merged 7 commits into from
Mar 21, 2024
Merged

Conversation

mari2895
Copy link
Collaborator

@mari2895 mari2895 commented Mar 20, 2024

PR Summary

PR Checklist

  • Adds a test for any bugs fixed. Adds tests for new features.
  • Format your changes by calling scripts/bash/format.sh.
  • Explain what you did.

@mari2895 mari2895 self-assigned this Mar 20, 2024
ReduceInGain<class internal_variables::GcovCool>(md, 0, 0);
};
auto Mdot400 = [](MeshData<Real> *md) {
Real rc = 0.04;
Copy link
Collaborator

Choose a reason for hiding this comment

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

Maybe make this a runtime set-able parameter. Might even be good to support an arbitrary number of these at different radii.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I am trying to imagine how is the second part going to work

Copy link
Collaborator

Choose a reason for hiding this comment

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

ParameterInput can read a std:::vector so you could do something like this:

auto mdot_radii = pin->GetOrAddVector("analysis", "mdot_radii", std::vector<Real>{0.04});
for (auto rc : mdot_radii) {
  auto mdot = [=](MeshData<Real> *md) { return History::CalculateMdot(md, rc, 0); };
  hst_vars.emplace_back(HistoryOutputVar(HstSum, mdot, "Mdot at r = " + std::to_string(rc)));
}

return ReduceInGain<class internal_variables::GcovHeat>(md, 0) -
ReduceInGain<class internal_variables::GcovCool>(md, 0);
return ReduceInGain<class internal_variables::GcovHeat>(md, 0, 0) -
ReduceInGain<class internal_variables::GcovCool>(md, 0, 0);
Copy link
Collaborator

Choose a reason for hiding this comment

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

Is the word class here actually required?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

deleted and still compiles and runs so is not required.

bool is_netheat = (v(b, iv::GcovHeat(), k, j, i) - v(b, iv::GcovCool(), k, j, i) >
1.e-8); // checks that in the gain region
auto analysis = pmb->packages.Get("analysis").get();
const Real outside_pns_threshold = analysis->Param<Real>("outside_pns_threshold");
Copy link
Collaborator

Choose a reason for hiding this comment

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

This needs to be moved outside of the kernel. It can't be called inside the par for.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

done

src/analysis/history.hpp Show resolved Hide resolved
src/analysis/history.cpp Outdated Show resolved Hide resolved
const F &f) {
Real div_mass_flux_integral;
div_mass_flux_integral =
-(f(b, 1, k, j, i + 1) - f(b, 1, k, j, i)) * coords.FaceArea<X1DIR>(k, j, i);
Copy link
Collaborator

Choose a reason for hiding this comment

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

Remove the minus signs here, and add them to the value in mdot as the divergence doesn't have a minus sign in its definition but mdot does, due to the orientation of surface normals.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

done

src/analysis/history.cpp Outdated Show resolved Hide resolved
Comment on lines 343 to 346
auto rad = pmb->packages.Get("radiation").get();
const parthenon::AllReduce<bool> *pdo_gain_reducer =
rad->MutableParam<parthenon::AllReduce<bool>>("do_gain_reducer");
const bool do_gain = pdo_gain_reducer->val;
Copy link
Collaborator

Choose a reason for hiding this comment

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

This can't be here---you can't access packages from within a par for on GPU. Move this above.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

done

Comment on lines 350 to 352
auto analysis = pmb->packages.Get("analysis").get();
const Real inside_pns_threshold =
analysis->Param<Real>("inside_pns_threshold");
Copy link
Collaborator

Choose a reason for hiding this comment

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

Can't be here. Move this out of the kernel.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

done

Comment on lines 9 to 15
Real outside_pns_threshold =
pin->GetOrAddReal("analysis", "outside_pns_threshold",
2.42e-5); // corresponds to entropy > 3 kb/baryon
Real inside_pns_threshold = pin->GetOrAddReal("analysis", "inside_pns_threshold",
0.008); // corresponds to r < 80 km
Real radius_cutoff_mdot =
pin->GetOrAddReal("analysis", "radius_cutoff_mdot", 0.04); // default 400km
Copy link
Collaborator

Choose a reason for hiding this comment

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

This stuff might better belong in the progenitor package.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

moved to progenitor pkg

ReduceInGain<class internal_variables::GcovCool>(md, 0, 0);
};
auto Mdot400 = [](MeshData<Real> *md) {
Real rc = 0.04;
Copy link
Collaborator

Choose a reason for hiding this comment

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

ParameterInput can read a std:::vector so you could do something like this:

auto mdot_radii = pin->GetOrAddVector("analysis", "mdot_radii", std::vector<Real>{0.04});
for (auto rc : mdot_radii) {
  auto mdot = [=](MeshData<Real> *md) { return History::CalculateMdot(md, rc, 0); };
  hst_vars.emplace_back(HistoryOutputVar(HstSum, mdot, "Mdot at r = " + std::to_string(rc)));
}

Copy link
Collaborator

@Yurlungur Yurlungur left a comment

Choose a reason for hiding this comment

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

This is looking good! A few more minor changes...

src/analysis/history.cpp Outdated Show resolved Hide resolved
src/analysis/history.cpp Outdated Show resolved Hide resolved
src/analysis/history.hpp Outdated Show resolved Hide resolved
src/analysis/history.hpp Show resolved Hide resolved
Comment on lines 94 to 97
Real LengthCGSToCode = units.GetLengthCGSToCode();
Real rc = 400e5 * LengthCGSToCode;
auto mdot_radii = pin->GetOrAddVector("progenitor", "mdot_radii",
std::vector<Real>{rc}); // default 400km
Copy link
Collaborator

Choose a reason for hiding this comment

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

I don't think you want to conversion here.

Suggested change
Real LengthCGSToCode = units.GetLengthCGSToCode();
Real rc = 400e5 * LengthCGSToCode;
auto mdot_radii = pin->GetOrAddVector("progenitor", "mdot_radii",
std::vector<Real>{rc}); // default 400km
Real LengthCGSToCode = units.GetLengthCGSToCode();
auto mdot_radii = pin->GetOrAddVector("progenitor", "mdot_radii",
std::vector<Real>{400}); // default 400km

Comment on lines 143 to 144
auto Mdot = [rc, LengthCGSToCode](MeshData<Real> *md) {
return History::CalculateMdot(md, rc, 0);
Copy link
Collaborator

Choose a reason for hiding this comment

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

Suggested change
auto Mdot = [rc, LengthCGSToCode](MeshData<Real> *md) {
return History::CalculateMdot(md, rc, 0);
auto rc_code = rc*1e5*LengthCGSToCode;
auto Mdot = [rc_code](MeshData<Real> *md) {
return History::CalculateMdot(md, rc_code, 0);

src/progenitor/progenitordata.cpp Outdated Show resolved Hide resolved
Copy link
Collaborator

@Yurlungur Yurlungur left a comment

Choose a reason for hiding this comment

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

As soon as tests pass here (and you confirm it works for you locally) I will merge.

@mari2895
Copy link
Collaborator Author

As soon as tests pass here (and you confirm it works for you locally) I will merge.

All tests have passed, compiles, runs, produces reasonable output.

@Yurlungur Yurlungur merged commit 799e02e into main Mar 21, 2024
2 checks passed
@Yurlungur Yurlungur deleted the mg/mdots branch March 21, 2024 16:23
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