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

Convert hdterm and hdterm3 to C++ #2625

Merged
merged 3 commits into from
Oct 10, 2023
Merged

Convert hdterm and hdterm3 to C++ #2625

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
35c669f
Convert hdterm and hdterm3 to C++
maxpkatz Oct 9, 2023
8271003
Merge branch 'development' into cxx-hdterm
maxpkatz Oct 9, 2023
de1e677
Fix typo
maxpkatz Oct 9, 2023
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
216 changes: 216 additions & 0 deletions Source/radiation/HABEC.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -367,6 +367,222 @@ namespace HABEC
}
});
}

AMREX_INLINE
void hdterm (Array4<Real> const dterm,
Array4<Real> const er,
const Box& reg,
int cdir, int bct, Real bcl,
Array4<Real const> const bcval,
Array4<int const> const mask,
Array4<Real const> const d,
const Real* dx)
{
Real h;

// Index shift for whichever edge we're checking.
// Negative means we're looking at the lo edge,
// positive means we're looking at the hi edge.
// The first group is for cell centers, the second
// group is for cell edges.

int icp = 0;
int jcp = 0;
int kcp = 0;

int iep = 0;
int jep = 0;
int kep = 0;

#if AMREX_SPACEDIM == 1
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 1) {
h = dx[0];
icp = 1;
iep = 1;
}
#elif AMREX_SPACEDIM == 2
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 2) {
h = dx[0];
icp = 1;
iep = 1;
}
else if (cdir == 1) {
h = dx[1];
jcp = -1;
}
else if (cdir == 3) {
h = dx[1];
jcp = 1;
jep = 1;
}
#else
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 3) {
h = dx[0];
icp = 1;
iep = 1;
}
else if (cdir == 1) {
h = dx[1];
jcp = -1;
}
else if (cdir == 4) {
h = dx[1];
jcp = 1;
jep = 1;
}
else if (cdir == 2) {
h = dx[2];
kcp = -1;
}
else if (cdir == 5) {
h = dx[2];
kcp = 1;
kep = 1;
}
#endif

amrex::LoopOnCpu(reg, [=] (int i, int j, int k) noexcept
{
if (mask.contains(i+icp,j+jcp,k+kcp)) {
if (mask(i+icp,j+jcp,k+kcp) > 0) {
if (bct == LO_DIRICHLET) {
Copy link
Member

Choose a reason for hiding this comment

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

do we want to throw an error if they come in with the wrong boundary type? That's the only real diff with the Fortran

Real d_sign = 1.0_rt;
if (iep != 0 || jep != 0 || kep != 0) {
// right edge
d_sign = -1.0_rt;
}
dterm(i+iep,j+jep,k+kep) = d(i+iep,j+jep,k+kep) * d_sign * (er(i,j,k) - bcval(i+icp,j+jcp,k+kcp)) / (0.5_rt * h + bcl);
}
}
}
});
}

AMREX_INLINE
void hdterm3 (Array4<Real> const dterm,
Array4<Real const> const er,
const Box& reg,
int cdir, int bctype,
Array4<int const> const tf,
Real bcl,
Array4<Real const> const bcval,
Array4<int const> const mask,
Array4<Real const> const d,
const Real* const dx)
{
Real h;

// Index shift for whichever edge we're checking.
// Negative means we're looking at the lo edge,
// positive means we're looking at the hi edge.
// The first group is for cell centers, the second
// group is for cell edges.

int icp = 0;
int jcp = 0;
Copy link
Member

Choose a reason for hiding this comment

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

we have this same code now in several places.. Maybe we should make an inline function with this and return things via structured binding? Doesn't need to be this PR

Copy link
Member Author

Choose a reason for hiding this comment

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

Yes, one of my plans after we finish the Fortran conversion is to do a bunch of cleanup on this since there's a bunch of similar-ish functions across the three Hypre ABEC implementations

int kcp = 0;

int iep = 0;
int jep = 0;
int kep = 0;

#if AMREX_SPACEDIM == 1
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 1) {
h = dx[0];
icp = 1;
iep = 1;
}
#elif AMREX_SPACEDIM == 2
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 2) {
h = dx[0];
icp = 1;
iep = 1;
}
else if (cdir == 1) {
h = dx[1];
jcp = -1;
}
else if (cdir == 3) {
h = dx[1];
jcp = 1;
jep = 1;
}
#else
if (cdir == 0) {
h = dx[0];
icp = -1;
}
else if (cdir == 3) {
h = dx[0];
icp = 1;
iep = 1;
}
else if (cdir == 1) {
h = dx[1];
jcp = -1;
}
else if (cdir == 4) {
h = dx[1];
jcp = 1;
jep = 1;
}
else if (cdir == 2) {
h = dx[2];
kcp = -1;
}
else if (cdir == 5) {
h = dx[2];
kcp = 1;
kep = 1;
}
#endif

amrex::LoopOnCpu(reg, [=] (int i, int j, int k) noexcept
{
if (mask.contains(i+icp,j+jcp,k+kcp)) {
if (mask(i+icp,j+jcp,k+kcp) > 0) {
int bct;
if (bctype == -1) {
bct = tf(i+icp,j+jcp,k+kcp);
}
else {
bct = bctype;
}
if (bct == LO_DIRICHLET) {
Real d_sign = 1.0_rt;
if (iep != 0 || jep != 0 || kep != 0) {
// right edge
d_sign = -1.0_rt;
}
dterm(i+iep,j+jep,k+kep) = d(i+iep,j+jep,k+kep) * d_sign * (er(i,j,k) - bcval(i+icp,j+jcp,k+kcp)) / (0.5_rt * h + bcl);
}
else if (bct == LO_NEUMANN && bcval(i+icp,j+jcp,k+kcp) == 0.0_rt) {
dterm(i+iep,j+jep,k+kep) = 0.0_rt;
}
}
}
});
}
} // namespace HABEC

#endif
133 changes: 0 additions & 133 deletions Source/radiation/HABEC_1D.F90
View file
Open in desktop

This file was deleted.

Loading