1#ifndef AMREX_ML_ABECLAPLACIAN_H_
2#define AMREX_ML_ABECLAPLACIAN_H_
3#include <AMReX_Config.H>
18 using FAB =
typename MF::fab_type;
19 using RT =
typename MF::value_type;
67 template <
typename T1,
typename T2,
68 std::enable_if_t<std::is_convertible_v<T1,typename MF::value_type> &&
69 std::is_convertible_v<T2,typename MF::value_type>,
82 template <
typename AMF,
83 std::enable_if_t<IsFabArray<AMF>::value &&
84 std::is_convertible_v<
typename AMF::value_type,
85 typename MF::value_type>,
87 void setACoeffs (
int amrlev,
const AMF& alpha);
99 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
112 template <
typename AMF,
113 std::enable_if_t<IsFabArray<AMF>::value &&
114 std::is_convertible_v<
typename AMF::value_type,
115 typename MF::value_type>,
117 void setBCoeffs (
int amrlev,
const Array<AMF const*,AMREX_SPACEDIM>& beta);
128 template <
typename T,
129 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
142 template <
typename T,
143 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
145 void setBCoeffs (
int amrlev, Vector<T>
const& beta);
157 void Fapply (
int amrlev,
int mglev, MF& out,
const MF& in)
const final;
158 void Fsmooth (
int amrlev,
int mglev, MF& sol,
const MF& rhs,
int redblack)
const final;
162 int face_only=0) const final;
164 void normalize (
int amrlev,
int mglev, MF& mf) const final;
168 [[nodiscard]] MF
const*
getACoeffs (
int amrlev,
int mglev)
const final
173 [[nodiscard]] std::unique_ptr<MLLinOpT<MF>>
makeNLinOp (
int )
const final;
188 static
void FFlux (
Box const& box, Real const* dxinv,
RT bscalar,
189 Array<
FAB const*, AMREX_SPACEDIM> const& bcoef,
190 Array<
FAB*,AMREX_SPACEDIM> const& flux,
191 FAB const& sol,
int face_only,
int ncomp);
218template <
typename MF>
226 define(a_geom, a_grids, a_dmap, a_info, a_factory, a_ncomp);
229template <
typename MF>
238 define(a_geom, a_grids, a_dmap, a_overset_mask, a_info, a_factory, a_ncomp);
243template <
typename MF>
253 this->m_ncomp = a_ncomp;
258template <
typename MF>
268 BL_PROFILE(
"MLABecLaplacian::define(overset)");
269 this->m_ncomp = a_ncomp;
274template <
typename MF>
278 m_a_coeffs.resize(this->m_num_amr_levels);
279 m_b_coeffs.resize(this->m_num_amr_levels);
280 for (
int amrlev = 0; amrlev < this->m_num_amr_levels; ++amrlev)
282 m_a_coeffs[amrlev].resize(this->m_num_mg_levels[amrlev]);
283 m_b_coeffs[amrlev].resize(this->m_num_mg_levels[amrlev]);
284 for (
int mglev = 0; mglev < this->m_num_mg_levels[amrlev]; ++mglev)
286 m_a_coeffs[amrlev][mglev].define
287 (this->m_grids[amrlev][mglev], this->m_dmap[amrlev][mglev],
288 1, 0,
MFInfo(), *(this->m_factory[amrlev][mglev]));
289 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
293 m_b_coeffs[amrlev][mglev][idim].define
294 (ba, this->m_dmap[amrlev][mglev], m_ncomp, 0,
MFInfo(),
295 *(this->m_factory[amrlev][mglev]));
301template <
typename MF>
302template <
typename T1,
typename T2,
303 std::enable_if_t<std::is_convertible_v<T1,typename MF::value_type> &&
304 std::is_convertible_v<T2,typename MF::value_type>,
int>>
310 if (m_a_scalar ==
RT(0.0)) {
311 for (
int amrlev = 0; amrlev < this->m_num_amr_levels; ++amrlev) {
312 m_a_coeffs[amrlev][0].setVal(
RT(0.0));
316 m_scalars_set =
true;
319template <
typename MF>
320template <
typename AMF,
321 std::enable_if_t<IsFabArray<AMF>::value &&
322 std::is_convertible_v<
typename AMF::value_type,
323 typename MF::value_type>,
int>>
328 "MLABecLaplacian::setACoeffs: alpha is supposed to be single component.");
329 m_a_coeffs[amrlev][0].LocalCopy(alpha, 0, 0, 1,
IntVect(0));
330 m_needs_update =
true;
334template <
typename MF>
336 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
int>>
340 m_a_coeffs[amrlev][0].setVal(
RT(alpha));
341 m_needs_update =
true;
346template <
typename MF>
347template <
typename AMF,
348 std::enable_if_t<IsFabArray<AMF>::value &&
349 std::is_convertible_v<
typename AMF::value_type,
350 typename MF::value_type>,
int>>
355 const int ncomp = this->getNComp();
357 if (beta[0]->
nComp() == ncomp) {
358 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
359 for (
int icomp = 0; icomp < ncomp; ++icomp) {
360 m_b_coeffs[amrlev][0][idim].LocalCopy(*beta[idim], icomp, icomp, 1,
IntVect(0));
364 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
365 for (
int icomp = 0; icomp < ncomp; ++icomp) {
366 m_b_coeffs[amrlev][0][idim].LocalCopy(*beta[idim], 0, icomp, 1,
IntVect(0));
370 m_needs_update =
true;
373template <
typename MF>
375 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
int>>
379 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
380 m_b_coeffs[amrlev][0][idim].setVal(
RT(beta));
382 m_needs_update =
true;
385template <
typename MF>
387 std::enable_if_t<std::is_convertible_v<T,typename MF::value_type>,
int>>
391 const int ncomp = this->getNComp();
392 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
393 for (
int icomp = 0; icomp < ncomp; ++icomp) {
394 m_b_coeffs[amrlev][0][idim].setVal(
RT(beta[icomp]));
397 m_needs_update =
true;
400template <
typename MF>
408#if (AMREX_SPACEDIM != 3)
409 applyMetricTermsCoeffs();
412 applyRobinBCTermsCoeffs();
416 update_singular_flags();
418 m_needs_update =
false;
421template <
typename MF>
425 BL_PROFILE(
"MLABecLaplacian::prepareForSolve()");
429#if (AMREX_SPACEDIM != 3)
430 applyMetricTermsCoeffs();
433 applyRobinBCTermsCoeffs();
437 update_singular_flags();
439 m_needs_update =
false;
442template <
typename MF>
446#if (AMREX_SPACEDIM != 3)
447 for (
int alev = 0; alev < this->m_num_amr_levels; ++alev)
450 this->applyMetricTerm(alev, mglev, m_a_coeffs[alev][mglev]);
451 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim)
453 this->applyMetricTerm(alev, mglev, m_b_coeffs[alev][mglev][idim]);
487template <
typename LP>
490 using RT =
typename LP::RT;
492 const int ncomp = linop.getNComp();
493 bool reset_alpha =
false;
494 if (linop.m_a_scalar == RT(0.0)) {
495 linop.m_a_scalar = RT(1.0);
498 const RT bovera = linop.m_b_scalar/linop.m_a_scalar;
502 "To reuse solver With Robin BC, one must re-call setScalars (and setACoeffs if the scalar is not zero)");
505 linop.m_scalars_set =
false;
506 linop.m_acoef_set =
false;
508 for (
int amrlev = 0; amrlev < linop.NAMRLevels(); ++amrlev) {
510 const Box& domain = linop.Geom(amrlev,mglev).Domain();
511 const RT dxi =
static_cast<RT
>(linop.Geom(amrlev,mglev).InvCellSize(0));
512 const RT dyi =
static_cast<RT
>((AMREX_SPACEDIM >= 2) ? linop.Geom(amrlev,mglev).InvCellSize(1) : Real(1.0));
513 const RT dzi =
static_cast<RT
>((AMREX_SPACEDIM == 3) ? linop.Geom(amrlev,mglev).InvCellSize(2) : Real(1.0));
516 linop.m_a_coeffs[amrlev][mglev].setVal(RT(0.0));
523#pragma omp parallel if (Gpu::notInLaunchRegion())
525 for (
MFIter mfi(linop.m_a_coeffs[amrlev][mglev], mfi_info); mfi.
isValid(); ++mfi)
527 const Box& vbx = mfi.validbox();
528 auto const& afab = linop.m_a_coeffs[amrlev][mglev].array(mfi);
529 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
530 auto const& bfab = linop.m_b_coeffs[amrlev][mglev][idim].const_array(mfi);
533 bool outside_domain_lo = !(domain.
contains(blo));
534 bool outside_domain_hi = !(domain.
contains(bhi));
535 if ((!outside_domain_lo) && (!outside_domain_hi)) {
continue; }
536 for (
int icomp = 0; icomp < ncomp; ++icomp) {
537 auto const& rbc = (*(linop.m_robin_bcval[amrlev]))[mfi].const_array(icomp*3);
538 if (linop.m_lobc_orig[icomp][idim] == LinOpBCType::Robin && outside_domain_lo)
541 RT fac = bovera*dxi*dxi;
544 RT B = (rbc(i,j,k,1)*dxi - rbc(i,j,k,0)*RT(0.5))
545 / (rbc(i,j,k,1)*dxi + rbc(i,j,k,0)*RT(0.5));
546 afab(i+1,j,k,icomp) += fac*bfab(i+1,j,k,icomp)*(RT(1.0)-B);
548 }
else if (idim == 1) {
549 RT fac = bovera*dyi*dyi;
552 RT B = (rbc(i,j,k,1)*dyi - rbc(i,j,k,0)*RT(0.5))
553 / (rbc(i,j,k,1)*dyi + rbc(i,j,k,0)*RT(0.5));
554 afab(i,j+1,k,icomp) += fac*bfab(i,j+1,k,icomp)*(RT(1.0)-B);
557 RT fac = bovera*dzi*dzi;
560 RT B = (rbc(i,j,k,1)*dzi - rbc(i,j,k,0)*RT(0.5))
561 / (rbc(i,j,k,1)*dzi + rbc(i,j,k,0)*RT(0.5));
562 afab(i,j,k+1,icomp) += fac*bfab(i,j,k+1,icomp)*(RT(1.0)-B);
566 if (linop.m_hibc_orig[icomp][idim] == LinOpBCType::Robin && outside_domain_hi)
569 RT fac = bovera*dxi*dxi;
572 RT B = (rbc(i,j,k,1)*dxi - rbc(i,j,k,0)*RT(0.5))
573 / (rbc(i,j,k,1)*dxi + rbc(i,j,k,0)*RT(0.5));
574 afab(i-1,j,k,icomp) += fac*bfab(i,j,k,icomp)*(RT(1.0)-B);
576 }
else if (idim == 1) {
577 RT fac = bovera*dyi*dyi;
580 RT B = (rbc(i,j,k,1)*dyi - rbc(i,j,k,0)*RT(0.5))
581 / (rbc(i,j,k,1)*dyi + rbc(i,j,k,0)*RT(0.5));
582 afab(i,j-1,k,icomp) += fac*bfab(i,j,k,icomp)*(RT(1.0)-B);
585 RT fac = bovera*dzi*dzi;
588 RT B = (rbc(i,j,k,1)*dzi - rbc(i,j,k,0)*RT(0.5))
589 / (rbc(i,j,k,1)*dzi + rbc(i,j,k,0)*RT(0.5));
590 afab(i,j,k-1,icomp) += fac*bfab(i,j,k,icomp)*(RT(1.0)-B);
601template <
typename MF>
605 if (this->hasRobinBC()) {
610template <
typename MF>
614 BL_PROFILE(
"MLABecLaplacian::averageDownCoeffs()");
616 for (
int amrlev = this->m_num_amr_levels-1; amrlev > 0; --amrlev)
618 auto& fine_a_coeffs = m_a_coeffs[amrlev];
619 auto& fine_b_coeffs = m_b_coeffs[amrlev];
621 averageDownCoeffsSameAmrLevel(amrlev, fine_a_coeffs, fine_b_coeffs);
622 averageDownCoeffsToCoarseAmrLevel(amrlev);
625 averageDownCoeffsSameAmrLevel(0, m_a_coeffs[0], m_b_coeffs[0]);
628template <
typename MF>
633 int nmglevs = a.
size();
634 for (
int mglev = 1; mglev < nmglevs; ++mglev)
636 IntVect ratio = (amrlev > 0) ?
IntVect(this->mg_coarsen_ratio) : this->mg_coarsen_ratio_vec[mglev-1];
638 if (m_a_scalar == 0.0)
640 a[mglev].setVal(
RT(0.0));
657 for (
int mglev = 1; mglev < nmglevs; ++mglev)
659 if (this->m_overset_mask[amrlev][mglev]) {
660 const RT fac =
static_cast<RT>(1 << mglev);
661 const RT osfac =
RT(2.0)*fac/(fac+
RT(1.0));
662 const int ncomp = this->getNComp();
664#pragma omp parallel if (Gpu::notInLaunchRegion())
669 Box const& ybx = mfi.nodaltilebox(1);,
670 Box const& zbx = mfi.nodaltilebox(2));
672 auto const& by =
b[mglev][1].array(mfi);,
673 auto const& bz =
b[mglev][2].array(mfi));
674 auto const& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
675#if defined(AMREX_USE_CUDA) && defined(_WIN32)
681#if (AMREX_SPACEDIM >= 2)
688#if (AMREX_SPACEDIM == 3)
715template <
typename MF>
719 auto& fine_a_coeffs = m_a_coeffs[flev ].back();
720 auto& fine_b_coeffs = m_b_coeffs[flev ].back();
721 auto& crse_a_coeffs = m_a_coeffs[flev-1].front();
722 auto& crse_b_coeffs = m_b_coeffs[flev-1].front();
724 if (m_a_scalar != 0.0) {
732 IntVect(this->mg_coarsen_ratio),
733 this->m_geom[flev-1][0]);
736template <
typename MF>
740 m_is_singular.clear();
741 m_is_singular.resize(this->m_num_amr_levels,
false);
742 auto itlo = std::find(this->m_lobc[0].
begin(), this->m_lobc[0].
end(), BCType::Dirichlet);
743 auto ithi = std::find(this->m_hibc[0].
begin(), this->m_hibc[0].
end(), BCType::Dirichlet);
744 if (itlo == this->m_lobc[0].
end() && ithi == this->m_hibc[0].
end())
746 for (
int alev = 0; alev < this->m_num_amr_levels; ++alev)
749 if (this->m_domain_covered[alev] && !this->m_overset_mask[alev][0])
751 if (m_a_scalar == Real(0.0))
753 m_is_singular[alev] =
true;
757 RT asum = m_a_coeffs[alev].back().sum(0,
IntVect(0));
758 RT amax = m_a_coeffs[alev].back().norminf(0,1,
IntVect(0));
759 m_is_singular[alev] = (std::abs(asum) <= amax *
RT(1.e-12));
765 if (!m_is_singular[0] && this->m_needs_coarse_data_for_bc &&
766 this->m_coarse_fine_bc_type == LinOpBCType::Neumann)
770 bool lev0_a_is_zero =
false;
771 if (m_a_scalar == Real(0.0)) {
772 lev0_a_is_zero =
true;
774 RT asum = m_a_coeffs[0].back().sum(0,
IntVect(0));
775 RT amax = m_a_coeffs[0].back().norminf(0,1,
IntVect(0));
776 bool a_is_almost_zero = std::abs(asum) <= amax *
RT(1.e-12);
777 if (a_is_almost_zero) { lev0_a_is_zero =
true; }
780 if (lev0_a_is_zero) {
781 auto bbox = this->m_grids[0][0].minimalBox();
782 for (
int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
783 if (this->m_lobc[0][idim] == LinOpBCType::Dirichlet) {
786 if (this->m_hibc[0][idim] == LinOpBCType::Dirichlet) {
790 if (this->m_geom[0][0].Domain().contains(bbox)) {
791 m_is_singular[0] =
true;
797template <
typename MF>
803 const MF& acoef = m_a_coeffs[amrlev][mglev];
804 AMREX_D_TERM(
const MF& bxcoef = m_b_coeffs[amrlev][mglev][0];,
805 const MF& bycoef = m_b_coeffs[amrlev][mglev][1];,
806 const MF& bzcoef = m_b_coeffs[amrlev][mglev][2];);
809 {
AMREX_D_DECL(
static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(0)),
810 static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(1)),
811 static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(2)))};
813 const RT ascalar = m_a_scalar;
814 const RT bscalar = m_b_scalar;
816 const int ncomp = this->getNComp();
820 const auto& xma = in.const_arrays();
821 const auto& yma = out.arrays();
822 const auto& ama = acoef.arrays();
824 const auto& byma = bycoef.const_arrays();,
825 const auto& bzma = bzcoef.const_arrays(););
826 if (this->m_overset_mask[amrlev][mglev]) {
827 const auto& osmma = this->m_overset_mask[amrlev][mglev]->const_arrays();
833 osmma[box_no], dxinv, ascalar, bscalar);
841 dxinv, ascalar, bscalar);
849#pragma omp parallel if (Gpu::notInLaunchRegion())
853 const Box& bx = mfi.tilebox();
854 const auto& xfab = in.array(mfi);
855 const auto& yfab = out.array(mfi);
856 const auto& afab = acoef.array(mfi);
858 const auto& byfab = bycoef.array(mfi);,
859 const auto& bzfab = bzcoef.array(mfi););
860 if (this->m_overset_mask[amrlev][mglev]) {
861 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
864 mlabeclap_adotx_os(i,j,k,n, yfab, xfab, afab,
AMREX_D_DECL(bxfab,byfab,bzfab),
865 osm, dxinv, ascalar, bscalar);
870 mlabeclap_adotx(i,j,k,n, yfab, xfab, afab,
AMREX_D_DECL(bxfab,byfab,bzfab),
871 dxinv, ascalar, bscalar);
878template <
typename MF>
884 bool regular_coarsening =
true;
885 if (amrlev == 0 && mglev > 0) {
886 regular_coarsening = this->mg_coarsen_ratio_vec[mglev-1] == this->mg_coarsen_ratio;
890 if (! this->m_use_gauss_seidel && regular_coarsening) {
891 Ax.define(sol.boxArray(), sol.DistributionMap(), sol.nComp(), 0);
892 Fapply(amrlev, mglev, Ax, sol);
895 const MF& acoef = m_a_coeffs[amrlev][mglev];
897 AMREX_D_TERM(
const MF& bxcoef = m_b_coeffs[amrlev][mglev][0];,
898 const MF& bycoef = m_b_coeffs[amrlev][mglev][1];,
899 const MF& bzcoef = m_b_coeffs[amrlev][mglev][2];);
900 const auto& undrrelxr = this->m_undrrelxr[amrlev][mglev];
901 const auto& maskvals = this->m_maskvals [amrlev][mglev];
905 const auto& f0 = undrrelxr[oitr()]; ++oitr;
906 const auto& f1 = undrrelxr[oitr()]; ++oitr;
907#if (AMREX_SPACEDIM > 1)
908 const auto& f2 = undrrelxr[oitr()]; ++oitr;
909 const auto& f3 = undrrelxr[oitr()]; ++oitr;
910#if (AMREX_SPACEDIM > 2)
911 const auto& f4 = undrrelxr[oitr()]; ++oitr;
912 const auto& f5 = undrrelxr[oitr()]; ++oitr;
918#if (AMREX_SPACEDIM > 1)
921#if (AMREX_SPACEDIM > 2)
927 const int nc = this->getNComp();
928 const Real* h = this->m_geom[amrlev][mglev].CellSize();
930 const RT dhy = m_b_scalar/
static_cast<RT>(h[1]*h[1]);,
931 const RT dhz = m_b_scalar/
static_cast<RT>(h[2]*h[2]));
932 const RT alpha = m_a_scalar;
936 && (this->m_overset_mask[amrlev][mglev] || regular_coarsening))
940#if (AMREX_SPACEDIM > 1)
943#if (AMREX_SPACEDIM > 2)
949 const auto& solnma = sol.arrays();
950 const auto& rhsma = rhs.const_arrays();
951 const auto& ama = acoef.const_arrays();
954 const auto& byma = bycoef.const_arrays();,
955 const auto& bzma = bzcoef.const_arrays(););
957 const auto& f0ma = f0.const_arrays();
958 const auto& f1ma = f1.const_arrays();
959#if (AMREX_SPACEDIM > 1)
960 const auto& f2ma = f2.const_arrays();
961 const auto& f3ma = f3.const_arrays();
962#if (AMREX_SPACEDIM > 2)
963 const auto& f4ma = f4.const_arrays();
964 const auto& f5ma = f5.const_arrays();
968 if (this->m_overset_mask[amrlev][mglev]) {
969 const auto& osmma = this->m_overset_mask[amrlev][mglev]->const_arrays();
970 if (this->m_use_gauss_seidel) {
974 Box vbx(ama[box_no]);
975 abec_gsrb_os(i,j,k,n, solnma[box_no], rhsma[box_no], alpha, ama[box_no],
982 osmma[box_no], vbx, redblack);
985 const auto& axma = Ax.const_arrays();
989 Box vbx(ama[box_no]);
990 abec_jacobi_os(i,j,k,n, solnma[box_no], rhsma[box_no], axma[box_no],
1001 }
else if (regular_coarsening) {
1002 if (this->m_use_gauss_seidel) {
1006 Box vbx(ama[box_no]);
1007 abec_gsrb(i,j,k,n, solnma[box_no], rhsma[box_no], alpha, ama[box_no],
1017 const auto& axma = Ax.const_arrays();
1021 Box vbx(ama[box_no]);
1022 abec_jacobi(i,j,k,n, solnma[box_no], rhsma[box_no], axma[box_no],
1042#pragma omp parallel if (Gpu::notInLaunchRegion())
1046 const auto& m0 = mm0.
array(mfi);
1047 const auto& m1 = mm1.
array(mfi);
1048#if (AMREX_SPACEDIM > 1)
1049 const auto& m2 = mm2.
array(mfi);
1050 const auto& m3 = mm3.
array(mfi);
1051#if (AMREX_SPACEDIM > 2)
1052 const auto& m4 = mm4.
array(mfi);
1053 const auto& m5 = mm5.
array(mfi);
1057 const Box& tbx = mfi.tilebox();
1058 const Box& vbx = mfi.validbox();
1059 const auto& solnfab = sol.array(mfi);
1060 const auto& rhsfab = rhs.const_array(mfi);
1061 const auto& afab = acoef.const_array(mfi);
1063 AMREX_D_TERM(
const auto& bxfab = bxcoef.const_array(mfi);,
1064 const auto& byfab = bycoef.const_array(mfi);,
1065 const auto& bzfab = bzcoef.const_array(mfi););
1067 const auto& f0fab = f0.const_array(mfi);
1068 const auto& f1fab = f1.const_array(mfi);
1069#if (AMREX_SPACEDIM > 1)
1070 const auto& f2fab = f2.const_array(mfi);
1071 const auto& f3fab = f3.const_array(mfi);
1072#if (AMREX_SPACEDIM > 2)
1073 const auto& f4fab = f4.const_array(mfi);
1074 const auto& f5fab = f5.const_array(mfi);
1078 if (this->m_overset_mask[amrlev][mglev]) {
1079 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
1080 if (this->m_use_gauss_seidel) {
1090 osm, vbx, redblack);
1093 const auto& axfab = Ax.const_array(mfi);
1107 }
else if (regular_coarsening) {
1108 if (this->m_use_gauss_seidel) {
1111 abec_gsrb(i,j,k,n, solnfab, rhsfab, alpha, afab,
1121 const auto& axfab = Ax.const_array(mfi);
1150template <
typename MF>
1158 const int mglev = 0;
1160 const Real* dxinv = this->m_geom[amrlev][mglev].InvCellSize();
1161 const int ncomp = this->getNComp();
1162 FFlux(box, dxinv, m_b_scalar,
1164 &(m_b_coeffs[amrlev][mglev][1][mfi]),
1165 &(m_b_coeffs[amrlev][mglev][2][mfi]))}},
1166 flux, sol, face_only, ncomp);
1169template <
typename MF>
1174 FAB const& sol,
int face_only,
int ncomp)
1177 const auto by = bcoef[1]->const_array();,
1178 const auto bz = bcoef[2]->const_array(););
1180 const auto& fyarr = flux[1]->array();,
1181 const auto& fzarr = flux[2]->array(););
1182 const auto& solarr = sol.array();
1186 RT fac = bscalar*
static_cast<RT>(dxinv[0]);
1188 int blen = box.
length(0);
1193#if (AMREX_SPACEDIM >= 2)
1194 fac = bscalar*
static_cast<RT>(dxinv[1]);
1202#if (AMREX_SPACEDIM == 3)
1203 fac = bscalar*
static_cast<RT>(dxinv[2]);
1214 RT fac = bscalar*
static_cast<RT>(dxinv[0]);
1220#if (AMREX_SPACEDIM >= 2)
1221 fac = bscalar*
static_cast<RT>(dxinv[1]);
1228#if (AMREX_SPACEDIM == 3)
1229 fac = bscalar*
static_cast<RT>(dxinv[2]);
1239template <
typename MF>
1245 const auto& acoef = m_a_coeffs[amrlev][mglev];
1246 AMREX_D_TERM(
const auto& bxcoef = m_b_coeffs[amrlev][mglev][0];,
1247 const auto& bycoef = m_b_coeffs[amrlev][mglev][1];,
1248 const auto& bzcoef = m_b_coeffs[amrlev][mglev][2];);
1251 {
AMREX_D_DECL(
static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(0)),
1252 static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(1)),
1253 static_cast<RT>(this->m_geom[amrlev][mglev].InvCellSize(2)))};
1255 const RT ascalar = m_a_scalar;
1256 const RT bscalar = m_b_scalar;
1258 const int ncomp = getNComp();
1262 const auto& ma = mf.arrays();
1263 const auto& ama = acoef.const_arrays();
1264 AMREX_D_TERM(
const auto& bxma = bxcoef.const_arrays();,
1265 const auto& byma = bycoef.const_arrays();,
1266 const auto& bzma = bzcoef.const_arrays(););
1272 dxinv, ascalar, bscalar);
1279#pragma omp parallel if (Gpu::notInLaunchRegion())
1283 const Box& bx = mfi.tilebox();
1284 const auto& fab = mf.array(mfi);
1285 const auto& afab = acoef.array(mfi);
1287 const auto& byfab = bycoef.array(mfi);,
1288 const auto& bzfab = bzcoef.array(mfi););
1293 dxinv, ascalar, bscalar);
1299template <
typename MF>
1303 bool support =
false;
1304 if (this->m_overset_mask[0][0]) {
1306 this->mg_domain_min_width)
1315template <
typename MF>
1316std::unique_ptr<MLLinOpT<MF>>
1319 if (this->m_overset_mask[0][0] ==
nullptr) {
return nullptr; }
1321 const Geometry& geom = this->m_geom[0].back();
1322 const BoxArray& ba = this->m_grids[0].back();
1325 std::unique_ptr<MLLinOpT<MF>>
r
1335 nop->setMaxOrder(this->maxorder);
1336 nop->setVerbose(this->
verbose);
1338 nop->setDomainBC(this->m_lobc, this->m_hibc);
1340 if (this->needsCoarseDataForBC())
1342 const Real* dx0 = this->m_geom[0][0].CellSize();
1343 RealVect fac(this->m_coarse_data_crse_ratio);
1346 nop->setCoarseFineBCLocation(cbloc);
1349 nop->setScalars(m_a_scalar, m_b_scalar);
1351 MF
const& alpha_bottom = m_a_coeffs[0].back();
1352 iMultiFab const& osm_bottom = *(this->m_overset_mask[0].back());
1353 const int ncomp = alpha_bottom.
nComp();
1354 MF alpha(ba, dm, ncomp, 0);
1356 RT a_max = alpha_bottom.norminf(0, ncomp,
IntVect(0),
true,
true);
1357 const int ncomp_b = m_b_coeffs[0].back()[0].nComp();
1359 RT by_max = m_b_coeffs[0].back()[1].
norminf(0,ncomp_b,
IntVect(0),
true,
true);,
1360 RT bz_max = m_b_coeffs[0].back()[2].
norminf(0,ncomp_b,
IntVect(0),
true,
true));
1365 RT huge_alpha =
RT(1.e30) *
1367 AMREX_D_DECL(std::abs(m_b_scalar)*bx_max*dxinv[0]*dxinv[0],
1368 std::abs(m_b_scalar)*by_max*dxinv[1]*dxinv[1],
1369 std::abs(m_b_scalar)*bz_max*dxinv[2]*dxinv[2]));
1374 auto const& ama = alpha.arrays();
1375 auto const& abotma = alpha_bottom.const_arrays();
1380 if (mma[box_no](i,j,k)) {
1381 ama[box_no](i,j,k,n) = abotma[box_no](i,j,k,n);
1383 ama[box_no](i,j,k,n) = huge_alpha;
1391#pragma omp parallel if (Gpu::notInLaunchRegion())
1394 Box const& bx = mfi.tilebox();
1395 auto const& a = alpha.array(mfi);
1396 auto const& abot = alpha_bottom.const_array(mfi);
1401 a(i,j,k,n) = abot(i,j,k,n);
1403 a(i,j,k,n) = huge_alpha;
1409 nop->setACoeffs(0, alpha);
1415template <
typename MF>
1419 if (this->m_overset_mask[0].back() ==
nullptr) {
return; }
1421 const int ncomp = dst.nComp();
1425 auto const& dstma = dst.arrays();
1426 auto const& srcma = src.const_arrays();
1427 auto const& mma = this->m_overset_mask[0].back()->const_arrays();
1431 if (mma[box_no](i,j,k)) {
1432 dstma[box_no](i,j,k,n) = srcma[box_no](i,j,k,n);
1434 dstma[box_no](i,j,k,n) =
RT(0.0);
1442#pragma omp parallel if (Gpu::notInLaunchRegion())
1445 Box const& bx = mfi.tilebox();
1446 auto const& dfab = dst.array(mfi);
1447 auto const& sfab = src.const_array(mfi);
1448 auto const& m = this->m_overset_mask[0].back()->const_array(mfi);
1452 dfab(i,j,k,n) = sfab(i,j,k,n);
1454 dfab(i,j,k,n) =
RT(0.0);
#define BL_PROFILE(a)
Definition AMReX_BLProfiler.H:551
#define AMREX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition AMReX_BLassert.H:49
#define AMREX_ASSERT_WITH_MESSAGE(EX, MSG)
Definition AMReX_BLassert.H:37
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition AMReX_BLassert.H:50
#define AMREX_GPU_LAUNCH_HOST_DEVICE_LAMBDA_RANGE(TN, TI, block)
Definition AMReX_GpuLaunchMacrosC.nolint.H:4
#define AMREX_HOST_DEVICE_FOR_3D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:106
#define AMREX_HOST_DEVICE_PARALLEL_FOR_4D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:111
#define AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM(a1, a2, a3, b1, b2, b3, c1, c2, c3)
Definition AMReX_GpuLaunch.nolint.H:29
#define AMREX_GPU_DEVICE
Definition AMReX_GpuQualifiers.H:18
Array4< Real > fine
Definition AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition AMReX_InterpFaceRegister.cpp:92
#define AMREX_LOOP_4D(bx, ncomp, i, j, k, n, block)
Definition AMReX_Loop.nolint.H:16
#define AMREX_D_TERM(a, b, c)
Definition AMReX_SPACE.H:129
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:550
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition AMReX_Box.H:146
AMREX_GPU_HOST_DEVICE bool contains(const IntVectND< dim > &p) const noexcept
Returns true if argument is contained within BoxND.
Definition AMReX_Box.H:204
const Real * InvCellSize() const noexcept
Returns the inverse cellsize for each coordinate direction.
Definition AMReX_CoordSys.H:82
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:41
int nComp() const noexcept
Return number of variables (aka components) associated with each point.
Definition AMReX_FabArrayBase.H:82
Array4< typename FabArray< FAB >::value_type const > const_array(const MFIter &mfi) const noexcept
Definition AMReX_FabArray.H:584
MultiArray4< typename FabArray< FAB >::value_type const > const_arrays() const noexcept
Definition AMReX_FabArray.H:646
Definition AMReX_FabFactory.H:50
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE constexpr IntVectND< dim > TheDimensionVector(int d) noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition AMReX_IntVect.H:691
Definition AMReX_MFIter.H:57
Box tilebox() const noexcept
Return the tile Box at the current index.
Definition AMReX_MFIter.cpp:385
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition AMReX_MFIter.H:141
Definition AMReX_MLABecLaplacian.H:15
void setACoeffs(int amrlev, const AMF &alpha)
Definition AMReX_MLABecLaplacian.H:325
RT getBScalar() const final
Definition AMReX_MLABecLaplacian.H:167
int getNComp() const override
Return number of components.
Definition AMReX_MLABecLaplacian.H:147
int m_ncomp
Definition AMReX_MLABecLaplacian.H:211
bool isSingular(int amrlev) const override
Is it singular on given AMR level?
Definition AMReX_MLABecLaplacian.H:155
LinOpBCType BCType
Definition AMReX_MLABecLaplacian.H:21
Array< MF const *, AMREX_SPACEDIM > getBCoeffs(int amrlev, int mglev) const final
Definition AMReX_MLABecLaplacian.H:170
void applyRobinBCTermsCoeffs()
Definition AMReX_MLABecLaplacian.H:603
void update_singular_flags()
Definition AMReX_MLABecLaplacian.H:738
typename MF::fab_type FAB
Definition AMReX_MLABecLaplacian.H:18
MLABecLaplacianT< MF > & operator=(const MLABecLaplacianT< MF > &)=delete
bool supportRobinBC() const noexcept override
Definition AMReX_MLABecLaplacian.H:207
RT m_b_scalar
Definition AMReX_MLABecLaplacian.H:194
typename MF::value_type RT
Definition AMReX_MLABecLaplacian.H:19
~MLABecLaplacianT() override
void prepareForSolve() override
Definition AMReX_MLABecLaplacian.H:423
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info=LPInfo(), const Vector< FabFactory< FAB > const * > &a_factory={}, int a_ncomp=1)
Definition AMReX_MLABecLaplacian.H:245
RT m_a_scalar
Definition AMReX_MLABecLaplacian.H:193
void define_ab_coeffs()
Definition AMReX_MLABecLaplacian.H:276
void averageDownCoeffsSameAmrLevel(int amrlev, Vector< MF > &a, Vector< Array< MF, AMREX_SPACEDIM > > &b)
Definition AMReX_MLABecLaplacian.H:630
void Fapply(int amrlev, int mglev, MF &out, const MF &in) const final
Definition AMReX_MLABecLaplacian.H:799
void averageDownCoeffsToCoarseAmrLevel(int flev)
Definition AMReX_MLABecLaplacian.H:717
bool m_scalars_set
Definition AMReX_MLABecLaplacian.H:198
Vector< int > m_is_singular
Definition AMReX_MLABecLaplacian.H:205
void setBCoeffs(int amrlev, const Array< AMF const *, AMREX_SPACEDIM > &beta)
Definition AMReX_MLABecLaplacian.H:352
void normalize(int amrlev, int mglev, MF &mf) const final
Divide mf by the diagonal component of the operator. Used by bicgstab.
Definition AMReX_MLABecLaplacian.H:1241
void averageDownCoeffs()
Definition AMReX_MLABecLaplacian.H:612
RT getAScalar() const final
Definition AMReX_MLABecLaplacian.H:166
bool m_needs_update
Definition AMReX_MLABecLaplacian.H:203
bool supportNSolve() const final
Definition AMReX_MLABecLaplacian.H:1301
std::unique_ptr< MLLinOpT< MF > > makeNLinOp(int) const final
Definition AMReX_MLABecLaplacian.H:1317
bool m_acoef_set
Definition AMReX_MLABecLaplacian.H:199
MLABecLaplacianT(MLABecLaplacianT< MF > &&)=delete
Vector< Vector< MF > > m_a_coeffs
Definition AMReX_MLABecLaplacian.H:195
void Fsmooth(int amrlev, int mglev, MF &sol, const MF &rhs, int redblack) const final
Definition AMReX_MLABecLaplacian.H:880
bool isBottomSingular() const override
Is the bottom of MG singular?
Definition AMReX_MLABecLaplacian.H:156
void applyMetricTermsCoeffs()
Definition AMReX_MLABecLaplacian.H:444
void copyNSolveSolution(MF &dst, MF const &src) const final
Definition AMReX_MLABecLaplacian.H:1417
MLABecLaplacianT()=default
typename MLLinOpT< MF >::Location Location
Definition AMReX_MLABecLaplacian.H:22
MF const * getACoeffs(int amrlev, int mglev) const final
Definition AMReX_MLABecLaplacian.H:168
void update() override
Update for reuse.
Definition AMReX_MLABecLaplacian.H:402
void FFlux(int amrlev, const MFIter &mfi, const Array< FAB *, AMREX_SPACEDIM > &flux, const FAB &sol, Location, int face_only=0) const final
Definition AMReX_MLABecLaplacian.H:1152
bool needsUpdate() const override
Does it need update if it's reused?
Definition AMReX_MLABecLaplacian.H:149
void setScalars(T1 a, T2 b) noexcept
Definition AMReX_MLABecLaplacian.H:306
MLABecLaplacianT(const MLABecLaplacianT< MF > &)=delete
Vector< Vector< Array< MF, AMREX_SPACEDIM > > > m_b_coeffs
Definition AMReX_MLABecLaplacian.H:196
Definition AMReX_MLCellABecLap.H:13
void define(const Vector< Geometry > &a_geom, const Vector< BoxArray > &a_grids, const Vector< DistributionMapping > &a_dmap, const LPInfo &a_info=LPInfo(), const Vector< FabFactory< FAB > const * > &a_factory={})
Definition AMReX_MLCellABecLap.H:94
void prepareForSolve() override
Definition AMReX_MLCellABecLap.H:247
void update() override
Update for reuse.
Definition AMReX_MLCellABecLap.H:240
static constexpr int mg_coarsen_ratio
Definition AMReX_MLLinOp.H:580
static constexpr int mg_box_min_width
Definition AMReX_MLLinOp.H:581
const MLLinOpT< MF > * m_parent
Definition AMReX_MLLinOp.H:596
Definition AMReX_MultiMask.H:18
MultiArray4< int const > const_arrays() const noexcept
Definition AMReX_MultiMask.H:48
Array4< int const > array(const MFIter &mfi) const noexcept
Definition AMReX_MultiMask.H:40
An Iterator over the Orientation of Faces of a Box.
Definition AMReX_Orientation.H:135
A Real vector in SpaceDim-dimensional space.
Definition AMReX_RealVect.H:32
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition AMReX_Vector.H:27
Long size() const noexcept
Definition AMReX_Vector.H:50
Definition AMReX_iMultiFab.H:32
void streamSynchronize() noexcept
Definition AMReX_GpuDevice.H:237
bool inLaunchRegion() noexcept
Definition AMReX_GpuControl.H:86
bool notInLaunchRegion() noexcept
Definition AMReX_GpuControl.H:87
void Max(KeyValuePair< K, V > &vi, MPI_Comm comm)
Definition AMReX_ParallelReduce.H:126
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition AMReX_ParallelContext.H:70
void applyRobinBCTermsCoeffs(LP &linop)
Definition AMReX_MLABecLaplacian.H:488
Definition AMReX_Amr.cpp:49
MF::value_type norminf(MF const &mf, int scomp, int ncomp, IntVect const &nghost, bool local=false)
Definition AMReX_FabArrayUtility.H:1883
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void overset_rescale_bcoef_x(Box const &box, Array4< T > const &bX, Array4< int const > const &osm, int ncomp, T osfac) noexcept
Definition AMReX_MLABecLap_1D_K.H:239
int nComp(FabArrayBase const &fa)
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition AMReX_CTOParallelForImpl.H:191
void average_down(const MultiFab &S_fine, MultiFab &S_crse, const Geometry &fgeom, const Geometry &cgeom, int scomp, int ncomp, int rr)
Definition AMReX_MultiFabUtil.cpp:314
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > surroundingNodes(const BoxND< dim > &b, int dir) noexcept
Returns a BoxND with NODE based coordinates in direction dir that encloses BoxND b....
Definition AMReX_Box.H:1399
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Returns a BoxND with different type.
Definition AMReX_Box.H:1435
void average_down_faces(const Vector< const MF * > &fine, const Vector< MF * > &crse, const IntVect &ratio, int ngcrse=0)
Average fine face-based FabArray onto crse face-based FabArray.
Definition AMReX_MultiFabUtil.H:851
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > adjCellLo(const BoxND< dim > &b, int dir, int len=1) noexcept
Returns the cell centered BoxND of length len adjacent to b on the low end along the coordinate direc...
Definition AMReX_Box.H:1591
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_adotx_os(int i, int, int, int n, Array4< T > const &y, Array4< T const > const &x, Array4< T const > const &a, Array4< T const > const &bX, Array4< int const > const &osm, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition AMReX_MLABecLap_1D_K.H:24
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_x(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, Array4< T const > const &bx, T fac, int ncomp) noexcept
Definition AMReX_MLABecLap_1D_K.H:56
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_z(Box const &box, Array4< T > const &fz, Array4< T const > const &sol, Array4< T const > const &bz, T fac, int ncomp) noexcept
Definition AMReX_MLABecLap_3D_K.H:163
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_y(Box const &box, Array4< T > const &fy, Array4< T const > const &sol, Array4< T const > const &by, T fac, int ncomp) noexcept
Definition AMReX_MLABecLap_2D_K.H:104
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition AMReX_Box.H:1890
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_gsrb_os(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Array4< int const > const &osm, Box const &vbox, int redblack) noexcept
Definition AMReX_MLABecLap_1D_K.H:121
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void overset_rescale_bcoef_y(Box const &box, Array4< T > const &bY, Array4< int const > const &osm, int ncomp, T osfac) noexcept
Definition AMReX_MLABecLap_2D_K.H:437
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_jacobi(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Box const &vbox) noexcept
Definition AMReX_MLABecLap_1D_K.H:160
IntVectND< AMREX_SPACEDIM > IntVect
Definition AMReX_BaseFwd.H:30
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > bdryLo(const BoxND< dim > &b, int dir, int len=1) noexcept
Returns the edge-centered BoxND (in direction dir) defining the low side of BoxND b.
Definition AMReX_Box.H:1502
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_gsrb(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Box const &vbox, int redblack) noexcept
Definition AMReX_MLABecLap_1D_K.H:87
AMREX_FORCE_INLINE void abec_gsrb_with_line_solve(Box const &, Array4< T > const &, Array4< T const > const &, T, Array4< T const > const &, T, Array4< T const > const &, Array4< int const > const &, Array4< int const > const &, Array4< T const > const &, Array4< T const > const &, Box const &, int, int) noexcept
Definition AMReX_MLABecLap_1D_K.H:223
std::array< T *, AMREX_SPACEDIM > GetArrOfPtrs(std::array< T, AMREX_SPACEDIM > &a) noexcept
Definition AMReX_Array.H:852
bool TilingIfNotGPU() noexcept
Definition AMReX_MFIter.H:12
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition AMReX_Box.H:1881
std::array< T const *, AMREX_SPACEDIM > GetArrOfConstPtrs(const std::array< T, AMREX_SPACEDIM > &a) noexcept
Definition AMReX_Array.H:864
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_yface(Box const &box, Array4< T > const &fy, Array4< T const > const &sol, Array4< T const > const &by, T fac, int ylen, int ncomp) noexcept
Definition AMReX_MLABecLap_2D_K.H:122
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_normalize(int i, int, int, int n, Array4< T > const &x, Array4< T const > const &a, Array4< T const > const &bX, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition AMReX_MLABecLap_1D_K.H:44
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr const T & max(const T &a, const T &b) noexcept
Definition AMReX_Algorithm.H:35
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_zface(Box const &box, Array4< T > const &fz, Array4< T const > const &sol, Array4< T const > const &bz, T fac, int zlen, int ncomp) noexcept
Definition AMReX_MLABecLap_3D_K.H:183
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > adjCellHi(const BoxND< dim > &b, int dir, int len=1) noexcept
Similar to adjCellLo but builds an adjacent BoxND on the high end.
Definition AMReX_Box.H:1612
int verbose
Definition AMReX_DistributionMapping.cpp:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void abec_jacobi_os(int i, int, int, int n, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T alpha, Array4< T const > const &a, T dhx, Array4< T const > const &bX, Array4< int const > const &m0, Array4< int const > const &m1, Array4< T const > const &f0, Array4< T const > const &f1, Array4< int const > const &osm, Box const &vbox) noexcept
Definition AMReX_MLABecLap_1D_K.H:189
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_flux_xface(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, Array4< T const > const &bx, T fac, int xlen, int ncomp) noexcept
Definition AMReX_MLABecLap_1D_K.H:72
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void overset_rescale_bcoef_z(Box const &box, Array4< T > const &bZ, Array4< int const > const &osm, int ncomp, T osfac) noexcept
Definition AMReX_MLABecLap_3D_K.H:777
std::array< T, N > Array
Definition AMReX_Array.H:24
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlabeclap_adotx(int i, int, int, int n, Array4< T > const &y, Array4< T const > const &x, Array4< T const > const &a, Array4< T const > const &bX, GpuArray< T, AMREX_SPACEDIM > const &dxinv, T alpha, T beta) noexcept
Definition AMReX_MLABecLap_1D_K.H:9
Definition AMReX_FabArrayCommI.H:896
Definition AMReX_Array.H:34
Definition AMReX_MLLinOp.H:35
Location
Definition AMReX_MLLinOp.H:87
FabArray memory allocation information.
Definition AMReX_FabArray.H:66
Definition AMReX_MFIter.H:20
MFItInfo & SetDynamic(bool f) noexcept
Definition AMReX_MFIter.H:34
MFItInfo & EnableTiling(const IntVect &ts=FabArrayBase::mfiter_tile_size) noexcept
Definition AMReX_MFIter.H:29