Skip to content

Commit

Permalink
Fix indexing error in computing the advective term for (AMReX-Fluids#104
Browse files Browse the repository at this point in the history 
)

non-conservative tracers. In practice, this makes little difference
because it's multiplied by div(Umac) which the MAC projection just
enforced to be numerically zero to within its tolerance.

Also make the density field non-constant and add a second,
non-conservative tracer in the tracer_adv_diff_cn regression test. This
will require creating new benchmarks for the 2D and 3D tests.

Co-authored-by: Ann Almgren <[email protected]>
  • Loading branch information
@cgilet @asalmgren
cgilet and asalmgren authored Jun 17, 2024
1 parent bc553d1 commit 0e232e8
Show file tree
Hide file tree
Showing 5 changed files with 55 additions and 23 deletions.
6 changes: 3 additions & 3 deletions src/convection/incflo_compute_advection_term.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -700,9 +700,9 @@ incflo::compute_convective_term (Vector<MultiFab*> const& conv_u,
// If convective, we define u dot grad trac = div (u trac) - trac div(u)
HydroUtils::ComputeConvectiveTerm(bx, m_ntrac, mfi,
tracer[lev]->array(mfi,0),
AMREX_D_DECL(face_x[lev].array(mfi),
face_y[lev].array(mfi),
face_z[lev].array(mfi)),
AMREX_D_DECL(face_x[lev].array(mfi,flux_comp),
face_y[lev].array(mfi,flux_comp),
face_z[lev].array(mfi,flux_comp)),
divu[lev].array(mfi),
update_arr,
get_tracer_iconserv_device_ptr(),
Expand Down
4 changes: 2 additions & 2 deletions src/prob/incflo_prob_I.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -84,14 +84,14 @@
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& problo,
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& probhi) const;

static void init_periodic_tracer (amrex::Box const& vbx, amrex::Box const& gbx,
void init_periodic_tracer (amrex::Box const& vbx, amrex::Box const& gbx,
amrex::Array4<amrex::Real> const& vel,
amrex::Array4<amrex::Real> const& density,
amrex::Array4<amrex::Real> const& tracer,
amrex::Box const& domain,
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& dx,
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& problo,
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& probhi);
amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& probhi) const;

void init_double_shear_layer (amrex::Box const& vbx, amrex::Box const& gbx,
amrex::Array4<amrex::Real> const& vel,
Expand Down
55 changes: 41 additions & 14 deletions src/prob/prob_init_fluid.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ void incflo::prob_init_fluid (int lev)
ld.tracer.array(mfi),
domain, dx, problo, probhi);
}
else if (12 == m_probtype)
else if (12 == m_probtype || 122 == m_probtype)
{
init_periodic_tracer(vbx, gbx,
ld.velocity.array(mfi),
Expand Down Expand Up @@ -809,30 +809,57 @@ void incflo::init_boussinesq_bubble (Box const& vbx, Box const& /*gbx*/,

void incflo::init_periodic_tracer (Box const& vbx, Box const& /*gbx*/,
Array4<Real> const& vel,
Array4<Real> const& /*density*/,
Array4<Real> const& density,
Array4<Real> const& tracer,
Box const& /*domain*/,
GpuArray<Real, AMREX_SPACEDIM> const& dx,
GpuArray<Real, AMREX_SPACEDIM> const& problo,
GpuArray<Real, AMREX_SPACEDIM> const& probhi)
GpuArray<Real, AMREX_SPACEDIM> const& probhi) const
{
Real L = probhi[0]-problo[0];
Real C = Real(2.0)*Real(3.1415926535897932) / L;
amrex::ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept

if (m_probtype == 12)
{
constexpr Real A = Real(1.0);
Real x = Real(i+0.5)*dx[0];
Real y = Real(j+0.5)*dx[1];
amrex::ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
constexpr Real A = Real(1.0);
Real x = Real(i+0.5)*dx[0];
Real y = Real(j+0.5)*dx[1];
#if (AMREX_SPACEDIM == 3)
Real z = Real(k+0.5)*dx[2];
Real z = Real(k+0.5)*dx[2];
#else
Real z = 0.0_rt;
Real z = 0.0_rt;
#endif
AMREX_D_TERM(vel(i,j,k,0) = Real(1.0);,
vel(i,j,k,1) = Real(0.1)*(std::sin(C*(x+z) - Real(0.00042)) + Real(1.0)) * std::exp(y);,
vel(i,j,k,2) = Real(0.1)*(std::sin(C*(x+y) - Real(0.00042)) + Real(1.0)) * std::exp(z););
tracer(i,j,k) = A *(std::sin(C*(y+z) - Real(0.00042)) + Real(1.0)) * std::exp(x);
});
AMREX_D_TERM(vel(i,j,k,0) = Real(1.0);,
vel(i,j,k,1) = Real(0.1)*(std::sin(C*(x+z) - Real(0.00042)) + Real(1.0)) * std::exp(y);,
vel(i,j,k,2) = Real(0.1)*(std::sin(C*(x+y) - Real(0.00042)) + Real(1.0)) * std::exp(z););
tracer(i,j,k,0) = A *(std::sin(C*(y+z) - Real(0.00042)) + Real(1.0)) * std::exp(x);
});
}
else if (m_probtype == 122)
{
Real B = Real(0.1);
amrex::ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
constexpr Real A = Real(1.0);
Real x = Real(i+0.5)*dx[0];
Real y = Real(j+0.5)*dx[1];
#if (AMREX_SPACEDIM == 3)
Real z = Real(k+0.5)*dx[2];
#else
Real z = 0.0_rt;
#endif
AMREX_D_TERM(vel(i,j,k,0) = Real(1.0);,
vel(i,j,k,1) = Real(0.1)*(std::sin(C*(x+z) - Real(0.00042)) + Real(1.0)) * std::exp(y);,
vel(i,j,k,2) = Real(0.1)*(std::sin(C*(x+y) - Real(0.00042)) + Real(1.0)) * std::exp(z););
density(i,j,k) = A + B*(x + y + z);
tracer(i,j,k,0) = A *(std::sin(C*(y+z) - Real(0.00042)) + Real(1.0)) * std::exp(x);
tracer(i,j,k,1) = A *(std::sin(C*(y+z) - Real(0.00042)) + Real(1.0)) * std::exp(x);
});
} else {
Abort("Unknow periodic tracer probtype");
}
}

void incflo::init_double_shear_layer (Box const& vbx, Box const& /*gbx*/,
Expand Down
6 changes: 4 additions & 2 deletions test_2d/benchmark.tracer_adv_diff_cn
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,9 @@ incflo.mu = 0.001 # Dynamic viscosity coefficient

incflo.constant_density = false #
incflo.advect_tracer = true #
incflo.mu_s = 0.1
incflo.ntrac = 2
incflo.mu_s = 0.1 0.1
incflo.trac_is_conservative = 1 0

incflo.diffusion_type = 1 # Crank-Nicolson (instead of Implicit)

Expand Down Expand Up @@ -71,7 +73,7 @@ incflo.write_eb_surface = false
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# INITIAL CONDITIONS #
#.......................................#
incflo.probtype = 12
incflo.probtype = 122
incflo.test_tracer_conservation = true

#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
Expand Down
7 changes: 5 additions & 2 deletions test_3d/benchmark.tracer_adv_diff_cn
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,10 @@ incflo.mu = 0.001 # Dynamic viscosity coefficient

incflo.constant_density = false #
incflo.advect_tracer = true #
incflo.mu_s = 0.1
incflo.ntrac = 2
incflo.mu_s = 0.1 0.1
incflo.trac_is_conservative = 1 0


incflo.diffusion_type = 1 # Crank-Nicolson (instead of Implicit)

Expand Down Expand Up @@ -75,7 +78,7 @@ incflo.write_eb_surface = false
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# INITIAL CONDITIONS #
#.......................................#
incflo.probtype = 12
incflo.probtype = 122
incflo.test_tracer_conservation = true

#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
Expand Down

0 comments on commit 0e232e8

Please sign in to comment.