Skip to content

Commit

Permalink
fix computation of avg_zmom (erf-model#1295)
Browse files Browse the repository at this point in the history 
* fix computation of avg_zmom

* fix white spaces

* add refluxing contribs from fast integrator
  • Loading branch information
@asalmgren
asalmgren authored Nov 14, 2023
1 parent 959cedb commit 0f73f17
Show file tree
Hide file tree
Showing 11 changed files with 291 additions and 141 deletions.
27 changes: 12 additions & 15 deletions Docs/sphinx_doc/MeshRefinement.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -124,7 +124,7 @@ computed by dividing the variable named rhotheta by the variable named density.
Coupling Types
--------------

ERF supports one-way, two-way, and "mixed" coupling; this is a run-time input
ERF supports one-way, two-way, and "mixed" coupling between levels; this is a run-time input

::

Expand All @@ -136,8 +136,9 @@ for the time advance of the fine solution . For cell-centered quantities,
and face-baced normal momenta on the coarse-fine interface, the coarse data is conservatively
interpolated to the fine mesh. The interpolated data is utilized to specify ghost cell data
(outside of the valid fine region) as well as specified and relaxation data inside the lateral boundaries
of the fine region. More specifically, a user may specify the total width of the interior
Dirichlet and relaxation region with ``erf.cf_width = <Int>`` (yellow + blue)
of the fine region. More specifically, similarly to how the lateral boundaries are treated,
a user may specify the total width of the interior Dirichlet and relaxation region with
``erf.cf_width = <Int>`` (yellow + blue)
and analogously the width of the interior Dirichlet region may be specified with
``erf.cf_set_width = <Int>`` (yellow).

Expand Down Expand Up @@ -173,7 +174,9 @@ where :math:`G` is the RHS of the evolution equations, :math:`\psi^{\prime}` is
relaxation, :math:`\psi^{FP}` is the fine data obtained from spatial and temporal interpolation of the
coarse data, and :math:`n` is the minimum number of grid points from a lateral boundary. The specified and
relaxation regions are applied to all dycore variables :math:`\left[\rho \; \rho\Theta \; U\; V\; W \right]`
on the fine mesh. Finally, we note that time dependent Dirichlet data, provided via an external boundary file,
on the fine mesh.

Finally, we note that time dependent Dirichlet data, provided via an external boundary file,
may be enforced on the lateral boundary conditions of the domain (coarsest mesh). For such cases,
the relaxation region width at the domain edges may be specified with ``erf.wrfbdy_width = <Int>``
(yellow + blue) while the interior Dirichlet region may be specified with ``erf.wrfbdy_set_width = <Int>``
Expand All @@ -187,17 +190,11 @@ the fine mesh communicates data back to the coarse mesh in two ways:

- A "reflux" operation is performed for all cell-centered data.

Because the normal momentum at the fine level is itself interpolated from the coarse level, the
difference between fluxes -- along the coarse-fine interfaces -- used to update the coarse data and fluxes
used to update the fine data is due to the difference in the averaging of the advected quantity to the face
where the flux is defined.

We note that both coupling schemes are conservative for mass because the fluxes for the continuity
equation are the momenta themselves, which are interpolated on faces at the coarse-fine interface. Other advected
quantities which are advanced in conservation form will lose conservation with one-way coupling.
Two-way coupling is conservative for these scalars as long as the refluxing operation is included with the
averaging down.

We define "mixed" coupling as using the two-way coupling algorithm for all cell-centered quantities except for
:math:`\rho` and :math:`\rho \theta.`

We note that all three coupling schemes are conservative for mass because the fluxes for the continuity
equation are the momenta themselves, which are interpolated on faces at the coarse-fine interface. Other advected
quantities which are advanced in conservation form will lose conservation with one-way coupling.
Two-way coupling ensures conservation of the advective contribution to all scalar updates but
does not account for loss of conservation due to diffusive or source terms.
20 changes: 10 additions & 10 deletions Source/Advection/AdvectionSrcForScalars.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,28 +72,28 @@ AdvectionSrcForScalarsVert(const amrex::Box& bx,
switch(vert_adv_type) {
case AdvType::Centered_2nd:
AdvectionSrcForScalarsWrapper<InterpType_H,CENTERED2>(bx, ncomp, icomp,
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Upwind_3rd:
AdvectionSrcForScalarsWrapper<InterpType_H,UPWIND3>(bx, ncomp, icomp,
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Centered_4th:
AdvectionSrcForScalarsWrapper<InterpType_H,CENTERED4>(bx, ncomp, icomp,
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Upwind_5th:
AdvectionSrcForScalarsWrapper<InterpType_H,UPWIND5>(bx, ncomp, icomp,
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Centered_6th:
AdvectionSrcForScalarsWrapper<InterpType_H,CENTERED6>(bx, ncomp, icomp,
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
break;
default:
AMREX_ASSERT_WITH_MESSAGE(false, "Unknown vertical advection scheme!");
Expand Down
20 changes: 10 additions & 10 deletions Source/Advection/AdvectionSrcForState.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -176,43 +176,43 @@ AdvectionSrcForScalars (const Box& bx, const int icomp, const int ncomp,
switch(horiz_adv_type) {
case AdvType::Centered_2nd:
AdvectionSrcForScalarsVert<CENTERED2>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
break;
case AdvType::Upwind_3rd:
AdvectionSrcForScalarsVert<UPWIND3>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
break;
case AdvType::Centered_4th:
AdvectionSrcForScalarsVert<CENTERED4>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
break;
case AdvType::Upwind_5th:
AdvectionSrcForScalarsVert<UPWIND5>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
break;
case AdvType::Centered_6th:
AdvectionSrcForScalarsVert<CENTERED6>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
avg_xmom, avg_ymom, avg_zmom, vert_adv_type);
break;
case AdvType::Weno_3:
AdvectionSrcForScalarsWrapper<WENO3,WENO3>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Weno_5:
AdvectionSrcForScalarsWrapper<WENO5,WENO5>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Weno_3Z:
AdvectionSrcForScalarsWrapper<WENO_Z3,WENO_Z3>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Weno_3MZQ:
AdvectionSrcForScalarsWrapper<WENO_MZQ3,WENO_MZQ3>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
avg_xmom, avg_ymom, avg_zmom);
break;
case AdvType::Weno_5Z:
AdvectionSrcForScalarsWrapper<WENO_Z5,WENO_Z5>(bx, ncomp, icomp, flx_arr, cell_prim,
avg_xmom, avg_ymom, avg_zmom);
avg_xmom, avg_ymom, avg_zmom);
break;
default:
AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
Expand Down
8 changes: 4 additions & 4 deletions Source/TimeIntegration/ERF_MRI.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,8 +18,8 @@ private:
std::function<void(T&, T&, T&, const amrex::Real, const amrex::Real, const amrex::Real, const int)> slow_rhs_pre;
std::function<void(T&, T&, T&, T&, const amrex::Real, const amrex::Real, const amrex::Real, const int)> slow_rhs_inc;
std::function<void(T&, T&, T&, T&, T&, const amrex::Real, const amrex::Real, const amrex::Real, const int )> slow_rhs_post;
std::function<void(int, int, T&, const T&, T&, T&, T&, const amrex::Real, const amrex::Real,
const amrex::Real, const amrex::Real)> fast_rhs;
std::function<void(int, int, int, T&, const T&, T&, T&, T&, const amrex::Real, const amrex::Real,
const amrex::Real, const amrex::Real)> fast_rhs;

/**
* \brief Integrator timestep size (Real)
Expand Down Expand Up @@ -133,7 +133,7 @@ public:
slow_rhs_post = F;
}

void set_fast_rhs (std::function<void(int, int, T&, const T&, T&, T&, T&,
void set_fast_rhs (std::function<void(int, int, int, T&, const T&, T&, T&, T&,
const amrex::Real, const amrex::Real,
const amrex::Real, const amrex::Real)> F)
{
Expand Down Expand Up @@ -281,7 +281,7 @@ public:
// *******************************************************************************
for (int ks = 0; ks < nsubsteps; ++ks)
{
fast_rhs(ks, nsubsteps, *F_slow, S_old, S_new, *S_sum, *S_scratch, dtau, inv_fac,
fast_rhs(ks, nsubsteps, nrk, *F_slow, S_old, S_new, *S_sum, *S_scratch, dtau, inv_fac,
time + ks*dtau, time + (ks+1) * dtau);

} // ks
Expand Down
Loading

0 comments on commit 0f73f17

Please sign in to comment.