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.

Sign up for GitHub

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

Avoid dividing by zero in constexpr linear algebra routines #1603

Merged
merged 2 commits into from
Jul 8, 2024
Merged

Avoid dividing by zero in constexpr linear algebra routines #1603

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
0380d0f
Avoid dividing by zero in constexpr linear algebra routines
yut23 Jul 5, 2024
b779e95
Have RHS::dgefa() return the info parameter directly
yut23 Jul 6, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 1 addition & 2 deletions integration/VODE/vode_dvjac.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -169,8 +169,7 @@ void dvjac (int& IERPJ, BurnT& state, DvodeT& vstate)
int IER{};

#ifdef NEW_NETWORK_IMPLEMENTATION
RHS::dgefa(vstate.jac);
IER = 0;
IER = RHS::dgefa(vstate.jac);
#else
if (integrator_rp::linalg_do_pivoting == 1) {
constexpr bool allow_pivot{true};
Expand Down
10 changes: 9 additions & 1 deletion networks/rhs.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -307,16 +307,22 @@ void dgesl (const RArray2D& a, RArray1D& b)
}

AMREX_GPU_HOST_DEVICE AMREX_INLINE
void dgefa (RArray2D& a)
int dgefa (RArray2D& a)
{

// LU factorization in-place without pivoting.

int info = 0;

amrex::constexpr_for<1, INT_NEQS>([&] (auto n1)
{
[[maybe_unused]] constexpr int k = n1;

// compute multipliers
if (a(k, k) == 0.0_rt) {
info = k;
return; // same as continue in a normal loop
}

amrex::Real t = -1.0_rt / a(k,k);
amrex::constexpr_for<k+1, INT_NEQS+1>([&] (auto n2)
Expand All @@ -340,6 +346,8 @@ void dgefa (RArray2D& a)
});
});
});

return info;
}

// Calculate the density dependence term for tabulated rates. The RHS has a term
Expand Down
Loading