Skip to content

Commit

Permalink
Call ca_shock from C++
Browse files Browse the repository at this point in the history
  • Loading branch information
maxpkatz committed Jun 29, 2019
1 parent 329e6d0 commit e3ba74e
Show file tree
Hide file tree
Showing 5 changed files with 49 additions and 103 deletions.
22 changes: 20 additions & 2 deletions Source/hydro/Castro_ctu_hydro.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -176,6 +176,26 @@ Castro::construct_ctu_hydro_source(Real time, Real dt)
shk.resize(obx, 1);
Elixir elix_shk = shk.elixir();
fab_size += shk.nBytes();

// Multidimensional shock detection
// Used for the hybrid Riemann solver

#ifdef SHOCK_VAR
bool compute_shock = true;
#else
bool compute_shock = false;
#endif

if (hybrid_riemann == 1 || compute_shock) {
#pragma gpu box(obx)
ca_shock(AMREX_INT_ANYD(obx.loVect()), AMREX_INT_ANYD(obx.hiVect()),
BL_TO_FORTRAN_ANYD(q[mfi]),
BL_TO_FORTRAN_ANYD(shk),
AMREX_REAL_ANYD(dx));
}
else {
shk.setVal(0.0);
}

qxm.resize(obx, NQ);
Elixir elix_qxm = qxm.elixir();
Expand Down Expand Up @@ -218,7 +238,6 @@ Castro::construct_ctu_hydro_source(Real time, Real dt)
BL_TO_FORTRAN_ANYD(flatn),
BL_TO_FORTRAN_ANYD(qaux[mfi]),
BL_TO_FORTRAN_ANYD(src_q[mfi]),
BL_TO_FORTRAN_ANYD(shk),
BL_TO_FORTRAN_ANYD(dq),
BL_TO_FORTRAN_ANYD(qxm),
BL_TO_FORTRAN_ANYD(qxp),
Expand Down Expand Up @@ -277,7 +296,6 @@ Castro::construct_ctu_hydro_source(Real time, Real dt)
BL_TO_FORTRAN_ANYD(flatn),
BL_TO_FORTRAN_ANYD(qaux[mfi]),
BL_TO_FORTRAN_ANYD(src_q[mfi]),
BL_TO_FORTRAN_ANYD(shk),
BL_TO_FORTRAN_ANYD(Ip),
BL_TO_FORTRAN_ANYD(Im),
BL_TO_FORTRAN_ANYD(Ip_src),
Expand Down
50 changes: 0 additions & 50 deletions Source/hydro/Castro_ctu_nd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,6 @@ subroutine ctu_ppm_states(lo, hi, &
flatn, f_lo, f_hi, &
qaux, qa_lo, qa_hi, &
srcQ, src_lo, src_hi, &
shk, sk_lo, sk_hi, &
Ip, Ip_lo, Ip_hi, &
Im, Im_lo, Im_hi, &
Ip_src, Ips_lo, Ips_hi, &
Expand Down Expand Up @@ -58,7 +57,6 @@ subroutine ctu_ppm_states(lo, hi, &
#else
use trace_ppm_module, only : trace_ppm
#endif
use advection_util_module, only : ca_shock
use prob_params_module, only : dg

implicit none
Expand All @@ -69,7 +67,6 @@ subroutine ctu_ppm_states(lo, hi, &
integer, intent(in) :: f_lo(3), f_hi(3)
integer, intent(in) :: qa_lo(3), qa_hi(3)
integer, intent(in) :: src_lo(3), src_hi(3)
integer, intent(in) :: sk_lo(3), sk_hi(3)
integer, intent(in) :: Ip_lo(3), Ip_hi(3)
integer, intent(in) :: Im_lo(3), Im_hi(3)
integer, intent(in) :: Ips_lo(3), Ips_hi(3)
Expand Down Expand Up @@ -100,7 +97,6 @@ subroutine ctu_ppm_states(lo, hi, &
real(rt), intent(in) :: flatn(f_lo(1):f_hi(1),f_lo(2):f_hi(2),f_lo(3):f_hi(3)) ! flattening parameter
real(rt), intent(in) :: srcQ(src_lo(1):src_hi(1),src_lo(2):src_hi(2),src_lo(3):src_hi(3),NQSRC) ! primitive variable source

real(rt), intent(inout) :: shk(sk_lo(1):sk_hi(1), sk_lo(2):sk_hi(2), sk_lo(3):sk_hi(3))
real(rt), intent(inout) :: Ip(Ip_lo(1):Ip_hi(1),Ip_lo(2):Ip_hi(2),Ip_lo(3):Ip_hi(3),1:3,NQ)
real(rt), intent(inout) :: Im(Im_lo(1):Im_hi(1),Im_lo(2):Im_hi(2),Im_lo(3):Im_hi(3),1:3,NQ)
real(rt), intent(inout) :: Ip_src(Ips_lo(1):Ips_hi(1),Ips_lo(2):Ips_hi(2),Ips_lo(3):Ips_hi(3),1:3,NQSRC)
Expand Down Expand Up @@ -130,31 +126,10 @@ subroutine ctu_ppm_states(lo, hi, &
logical :: source_nonzero(NQSRC)
logical :: reconstruct_state(NQ)

logical :: compute_shock

!$gpu

hdt = HALF*dt

! multidimensional shock detection

#ifdef SHOCK_VAR
compute_shock = .true.
#else
compute_shock = .false.
#endif

! multidimensional shock detection -- this will be used to do the
! hybrid Riemann solver
if (hybrid_riemann == 1 .or. compute_shock) then
call ca_shock(lo, hi, &
q, qd_lo, qd_hi, &
shk, sk_lo, sk_hi, &
dx)
else
shk(lo(1):hi(1),lo(2):hi(2),lo(3):hi(3)) = ZERO
endif

! we don't need to reconstruct all of the NQ state variables,
! depending on how we are tracing
reconstruct_state(:) = .true.
Expand Down Expand Up @@ -385,7 +360,6 @@ subroutine ctu_plm_states(lo, hi, &
flatn, f_lo, f_hi, &
qaux, qa_lo, qa_hi, &
srcQ, src_lo, src_hi, &
shk, sk_lo, sk_hi, &
dq, dq_lo, dq_hi, &
qxm, qxm_lo, qxm_hi, &
qxp, qxp_lo, qxp_hi, &
Expand Down Expand Up @@ -418,7 +392,6 @@ subroutine ctu_plm_states(lo, hi, &
plm_iorder, use_pslope, hybrid_riemann
use trace_plm_module, only : trace_plm
use slope_module, only : uslope, pslope
use advection_util_module, only : ca_shock
use prob_params_module, only : dg

implicit none
Expand All @@ -429,7 +402,6 @@ subroutine ctu_plm_states(lo, hi, &
integer, intent(in) :: f_lo(3), f_hi(3)
integer, intent(in) :: qa_lo(3), qa_hi(3)
integer, intent(in) :: src_lo(3), src_hi(3)
integer, intent(in) :: sk_lo(3), sk_hi(3)
integer, intent(in) :: dq_lo(3), dq_hi(3)
integer, intent(in) :: qxm_lo(3), qxm_hi(3)
integer, intent(in) :: qxp_lo(3), qxp_hi(3)
Expand All @@ -453,7 +425,6 @@ subroutine ctu_plm_states(lo, hi, &
real(rt), intent(in) :: flatn(f_lo(1):f_hi(1),f_lo(2):f_hi(2),f_lo(3):f_hi(3)) ! flattening parameter
real(rt), intent(in) :: srcQ(src_lo(1):src_hi(1),src_lo(2):src_hi(2),src_lo(3):src_hi(3),NQSRC) ! primitive variable source

real(rt), intent(inout) :: shk(sk_lo(1):sk_hi(1), sk_lo(2):sk_hi(2), sk_lo(3):sk_hi(3))
real(rt), intent(inout) :: dq(dq_lo(1):dq_hi(1), dq_lo(2):dq_hi(2), dq_lo(3):dq_hi(3), NQ)

real(rt), intent(inout) :: qxm(qxm_lo(1):qxm_hi(1), qxm_lo(2):qxm_hi(2), qxm_lo(3):qxm_hi(3), NQ)
Expand All @@ -474,31 +445,10 @@ subroutine ctu_plm_states(lo, hi, &

logical :: reconstruct_state(NQ)

logical :: compute_shock

!$gpu

hdt = HALF*dt

! multidimensional shock detection

#ifdef SHOCK_VAR
compute_shock = .true.
#else
compute_shock = .false.
#endif

! multidimensional shock detection -- this will be used to do the
! hybrid Riemann solver
if (hybrid_riemann == 1 .or. compute_shock) then
call ca_shock(lo, hi, &
q, qd_lo, qd_hi, &
shk, sk_lo, sk_hi, &
dx)
else
shk(lo(1):hi(1),lo(2):hi(2),lo(3):hi(3)) = ZERO
endif

! we don't need to reconstruct all of the NQ state variables,
! depending on how we are tracing
reconstruct_state(:) = .true.
Expand Down
10 changes: 6 additions & 4 deletions Source/hydro/Castro_hydro_F.H
Original file line number Diff line number Diff line change
Expand Up @@ -162,6 +162,12 @@ extern "C"
(const int* lo, const int* hi,
BL_FORT_FAB_ARG_3D(flux));

void ca_shock
(const int* lo, const int* hi,
const BL_FORT_FAB_ARG_3D(q),
BL_FORT_FAB_ARG_3D(shk),
const amrex::Real* dx);

void ca_ctu_update
(const int* lo, const int* hi,
const int* is_finest_level,
Expand Down Expand Up @@ -225,7 +231,6 @@ extern "C"
const BL_FORT_FAB_ARG_3D(flatn),
const BL_FORT_FAB_ARG_3D(qaux),
const BL_FORT_FAB_ARG_3D(srcQ),
const BL_FORT_FAB_ARG_3D(shk),
BL_FORT_FAB_ARG_3D(Ip),
BL_FORT_FAB_ARG_3D(Im),
BL_FORT_FAB_ARG_3D(Ip_src),
Expand Down Expand Up @@ -257,7 +262,6 @@ extern "C"
const BL_FORT_FAB_ARG_3D(flatn),
const BL_FORT_FAB_ARG_3D(qaux),
const BL_FORT_FAB_ARG_3D(srcQ),
const BL_FORT_FAB_ARG_3D(shk),
BL_FORT_FAB_ARG_3D(dq),
BL_FORT_FAB_ARG_3D(qxm),
BL_FORT_FAB_ARG_3D(qxp),
Expand Down Expand Up @@ -415,7 +419,6 @@ extern "C"
(const int* lo, const int* hi,
BL_FORT_FAB_ARG_3D(q),
const BL_FORT_FAB_ARG_3D(flatn),
BL_FORT_FAB_ARG_3D(shk),
BL_FORT_FAB_ARG_3D(dq),
BL_FORT_FAB_ARG_3D(qm),
BL_FORT_FAB_ARG_3D(qp),
Expand All @@ -425,7 +428,6 @@ extern "C"
(const int* lo, const int* hi,
BL_FORT_FAB_ARG_3D(q),
const BL_FORT_FAB_ARG_3D(flatn),
BL_FORT_FAB_ARG_3D(shk),
BL_FORT_FAB_ARG_3D(qm),
BL_FORT_FAB_ARG_3D(qp),
const amrex::Real* dx);
Expand Down
22 changes: 20 additions & 2 deletions Source/hydro/Castro_mol_hydro.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -180,6 +180,26 @@ Castro::construct_mol_hydro_source(Real time, Real dt, MultiFab& A_update)

shk.resize(obx, 1);
Elixir elix_shk = shk.elixir();

// Multidimensional shock detection
// Used for the hybrid Riemann solver

#ifdef SHOCK_VAR
bool compute_shock = true;
#else
bool compute_shock = false;
#endif

if (hybrid_riemann == 1 || compute_shock) {
#pragma gpu box(obx)
ca_shock(AMREX_INT_ANYD(obx.loVect()), AMREX_INT_ANYD(obx.hiVect()),
BL_TO_FORTRAN_ANYD(q[mfi]),
BL_TO_FORTRAN_ANYD(shk),
AMREX_REAL_ANYD(dx));
}
else {
shk.setVal(0.0);
}

qm.resize(tbx, NQ*AMREX_SPACEDIM);
Elixir elix_qm = qm.elixir();
Expand All @@ -199,7 +219,6 @@ Castro::construct_mol_hydro_source(Real time, Real dt, MultiFab& A_update)
(AMREX_INT_ANYD(obx.loVect()), AMREX_INT_ANYD(obx.hiVect()),
BL_TO_FORTRAN_ANYD(q[mfi]),
BL_TO_FORTRAN_ANYD(flatn),
BL_TO_FORTRAN_ANYD(shk),
BL_TO_FORTRAN_ANYD(dq),
BL_TO_FORTRAN_ANYD(qm),
BL_TO_FORTRAN_ANYD(qp),
Expand All @@ -212,7 +231,6 @@ Castro::construct_mol_hydro_source(Real time, Real dt, MultiFab& A_update)
(AMREX_INT_ANYD(obx.loVect()), AMREX_INT_ANYD(obx.hiVect()),
BL_TO_FORTRAN_ANYD(q[mfi]),
BL_TO_FORTRAN_ANYD(flatn),
BL_TO_FORTRAN_ANYD(shk),
BL_TO_FORTRAN_ANYD(qm),
BL_TO_FORTRAN_ANYD(qp),
AMREX_REAL_ANYD(dx));
Expand Down
48 changes: 3 additions & 45 deletions Source/hydro/Castro_mol_nd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -4,15 +4,14 @@
subroutine ca_mol_plm_reconstruct(lo, hi, &
q, q_lo, q_hi, &
flatn, fl_lo, fl_hi, &
shk, shk_lo, shk_hi, &
dq, dq_lo, dq_hi, &
qm, qm_lo, qm_hi, &
qp, qp_lo, qp_hi, &
dx) bind(C, name="ca_mol_plm_reconstruct")

use amrex_error_module
use meth_params_module, only : NQ, NVAR, NGDNV, GDPRES, &
UTEMP, USHK, UMX, &
UTEMP, UMX, &
use_flattening, QPRES, &
QTEMP, QFS, QFX, QREINT, QRHO, &
first_order_hydro, hybrid_riemann, &
Expand All @@ -24,48 +23,28 @@ subroutine ca_mol_plm_reconstruct(lo, hi, &
use eos_module, only : eos
use network, only : nspec, naux
use prob_params_module, only : dg, coord_type
use advection_util_module, only : ca_shock

implicit none

integer, intent(in) :: lo(3), hi(3)
integer, intent(in) :: q_lo(3), q_hi(3)
integer, intent(in) :: fl_lo(3), fl_hi(3)
integer, intent(in) :: shk_lo(3), shk_hi(3)
integer, intent(in) :: dq_lo(3), dq_hi(3)
integer, intent(in) :: qm_lo(3), qm_hi(3)
integer, intent(in) :: qp_lo(3), qp_hi(3)

real(rt), intent(inout) :: q(q_lo(1):q_hi(1), q_lo(2):q_hi(2), q_lo(3):q_hi(3), NQ)
real(rt), intent(in) :: flatn(fl_lo(1):fl_hi(1), fl_lo(2):fl_hi(2), fl_lo(3):fl_hi(3))
real(rt), intent(inout) :: shk(shk_lo(1):shk_hi(1), shk_lo(2):shk_hi(2), shk_lo(3):shk_hi(3))
real(rt), intent(inout) :: dq(dq_lo(1):dq_hi(1), dq_lo(2):dq_hi(2), dq_lo(3):dq_hi(3), NQ)
real(rt), intent(inout) :: qm(qm_lo(1):qm_hi(1), qm_lo(2):qm_hi(2), qm_lo(3):qm_hi(3), NQ, AMREX_SPACEDIM)
real(rt), intent(inout) :: qp(qp_lo(1):qp_hi(1), qp_lo(2):qp_hi(2), qp_lo(3):qp_hi(3), NQ, AMREX_SPACEDIM)
real(rt), intent(in) :: dx(3)


integer :: idir, i, j, k, n
logical :: compute_shock
type (eos_t) :: eos_state

!$gpu

#ifdef SHOCK_VAR
compute_shock = .true.
#else
compute_shock = .false.
#endif

if (hybrid_riemann == 1 .or. compute_shock) then
call ca_shock(lo, hi, &
q, q_lo, q_hi, &
shk, shk_lo, shk_hi, &
dx)
else
shk(lo(1):hi(1),lo(2):hi(2),lo(3):hi(3)) = ZERO
endif

do idir = 1, AMREX_SPACEDIM

do n = 1, NQ
Expand Down Expand Up @@ -175,14 +154,13 @@ end subroutine ca_mol_plm_reconstruct
subroutine ca_mol_ppm_reconstruct(lo, hi, &
q, q_lo, q_hi, &
flatn, fl_lo, fl_hi, &
shk, shk_lo, shk_hi, &
qm, qm_lo, qm_hi, &
qp, qp_lo, qp_hi, &
dx) bind(C, name="ca_mol_ppm_reconstruct")

use amrex_error_module
use meth_params_module, only : NQ, NVAR, NGDNV, GDPRES, &
UTEMP, USHK, UMX, &
UTEMP, UMX, &
use_flattening, QPRES, &
QTEMP, QFS, QFX, QREINT, QRHO, &
first_order_hydro, hybrid_riemann, &
Expand All @@ -194,46 +172,26 @@ subroutine ca_mol_ppm_reconstruct(lo, hi, &
use eos_module, only : eos
use network, only : nspec, naux
use prob_params_module, only : dg, coord_type
use advection_util_module, only : ca_shock

implicit none

integer, intent(in) :: lo(3), hi(3)
integer, intent(in) :: q_lo(3), q_hi(3)
integer, intent(in) :: fl_lo(3), fl_hi(3)
integer, intent(in) :: shk_lo(3), shk_hi(3)
integer, intent(in) :: qm_lo(3), qm_hi(3)
integer, intent(in) :: qp_lo(3), qp_hi(3)

real(rt), intent(inout) :: q(q_lo(1):q_hi(1), q_lo(2):q_hi(2), q_lo(3):q_hi(3), NQ)
real(rt), intent(in) :: flatn(fl_lo(1):fl_hi(1), fl_lo(2):fl_hi(2), fl_lo(3):fl_hi(3))
real(rt), intent(inout) :: shk(shk_lo(1):shk_hi(1), shk_lo(2):shk_hi(2), shk_lo(3):shk_hi(3))
real(rt), intent(inout) :: qm(qm_lo(1):qm_hi(1), qm_lo(2):qm_hi(2), qm_lo(3):qm_hi(3), NQ, AMREX_SPACEDIM)
real(rt), intent(inout) :: qp(qp_lo(1):qp_hi(1), qp_lo(2):qp_hi(2), qp_lo(3):qp_hi(3), NQ, AMREX_SPACEDIM)
real(rt), intent(in) :: dx(3)

integer :: idir, i, j, k, n
logical :: compute_shock
type (eos_t) :: eos_state

!$gpu

#ifdef SHOCK_VAR
compute_shock = .true.
#else
compute_shock = .false.
#endif

if (hybrid_riemann == 1 .or. compute_shock) then
call ca_shock(lo, hi, &
q, q_lo, q_hi, &
shk, shk_lo, shk_hi, &
dx)
else
shk(lo(1):hi(1),lo(2):hi(2),lo(3):hi(3)) = ZERO
endif


do idir = 1, AMREX_SPACEDIM
call ca_ppm_reconstruct(lo, hi, 1, idir, &
q, q_lo, q_hi, NQ, 1, NQ, &
Expand Down Expand Up @@ -310,7 +268,7 @@ subroutine ca_mol_consup(lo, hi, &

use amrex_error_module
use meth_params_module, only : NQ, NVAR, NGDNV, GDPRES, &
UTEMP, USHK, UMX, &
UTEMP, UMX, &
use_flattening, QPRES, &
QTEMP, QFS, QFX, QREINT, QRHO, &
first_order_hydro, difmag, hybrid_riemann, &
Expand Down

0 comments on commit e3ba74e

Please sign in to comment.