Block-Structured AMR Software Framework
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
AMReX_MLPoisson.H
Go to the documentation of this file.
1#ifndef AMREX_MLPOISSON_H_
2#define AMREX_MLPOISSON_H_
3#include <AMReX_Config.H>
4
6#include <AMReX_MLPoisson_K.H>
8
9namespace amrex {
10
11// del dot grad phi
12
13template <typename MF>
15 : public MLCellABecLapT<MF>
16{
17public:
18
19 using FAB = typename MF::fab_type;
20 using RT = typename MF::value_type;
21
22 using BCType = LinOpBCType;
24
25 MLPoissonT () = default;
26 MLPoissonT (const Vector<Geometry>& a_geom,
27 const Vector<BoxArray>& a_grids,
28 const Vector<DistributionMapping>& a_dmap,
29 const LPInfo& a_info = LPInfo(),
30 const Vector<FabFactory<FAB> const*>& a_factory = {});
31 MLPoissonT (const Vector<Geometry>& a_geom,
32 const Vector<BoxArray>& a_grids,
33 const Vector<DistributionMapping>& a_dmap,
34 const Vector<iMultiFab const*>& a_overset_mask, // 1: unknown, 0: known
35 const LPInfo& a_info = LPInfo(),
36 const Vector<FabFactory<FAB> const*>& a_factory = {});
37 ~MLPoissonT () override;
38
39 MLPoissonT (const MLPoissonT<MF>&) = delete;
43
44 void define (const Vector<Geometry>& a_geom,
45 const Vector<BoxArray>& a_grids,
46 const Vector<DistributionMapping>& a_dmap,
47 const LPInfo& a_info = LPInfo(),
48 const Vector<FabFactory<FAB> const*>& a_factory = {});
49
50 void define (const Vector<Geometry>& a_geom,
51 const Vector<BoxArray>& a_grids,
52 const Vector<DistributionMapping>& a_dmap,
53 const Vector<iMultiFab const*>& a_overset_mask,
54 const LPInfo& a_info = LPInfo(),
55 const Vector<FabFactory<FAB> const*>& a_factory = {});
56
57 void prepareForSolve () final;
58 [[nodiscard]] bool isSingular (int amrlev) const final { return m_is_singular[amrlev]; }
59 [[nodiscard]] bool isBottomSingular () const final { return m_is_singular[0]; }
60 void Fapply (int amrlev, int mglev, MF& out, const MF& in) const final;
61 void Fsmooth (int amrlev, int mglev, MF& sol, const MF& rhs, int redblack) const final;
62 void FFlux (int amrlev, const MFIter& mfi,
64 const FAB& sol, Location loc, int face_only=0) const final;
65
66 void normalize (int amrlev, int mglev, MF& mf) const final;
67
68 [[nodiscard]] RT getAScalar () const final { return RT(0.0); }
69 [[nodiscard]] RT getBScalar () const final { return RT(-1.0); }
70 [[nodiscard]] MF const* getACoeffs (int /*amrlev*/, int /*mglev*/) const final { return nullptr; }
71 [[nodiscard]] Array<MF const*,AMREX_SPACEDIM> getBCoeffs (int /*amrlev*/, int /*mglev*/) const final
72 { return {{ AMREX_D_DECL(nullptr,nullptr,nullptr)}}; }
73
74 [[nodiscard]] std::unique_ptr<MLLinOpT<MF>> makeNLinOp (int grid_size) const final;
75
76 [[nodiscard]] bool supportNSolve () const final;
77
78 void copyNSolveSolution (MF& dst, MF const& src) const final;
79
81 void get_dpdn_on_domain_faces (Array<MF*,AMREX_SPACEDIM> const& dpdn,
82 MF const& phi);
83
84private:
85
87};
88
89template <typename MF>
90MLPoissonT<MF>::MLPoissonT (const Vector<Geometry>& a_geom,
91 const Vector<BoxArray>& a_grids,
92 const Vector<DistributionMapping>& a_dmap,
93 const LPInfo& a_info,
94 const Vector<FabFactory<FAB> const*>& a_factory)
95{
96 define(a_geom, a_grids, a_dmap, a_info, a_factory);
97}
98
99template <typename MF>
101 const Vector<BoxArray>& a_grids,
102 const Vector<DistributionMapping>& a_dmap,
103 const Vector<iMultiFab const*>& a_overset_mask,
104 const LPInfo& a_info,
105 const Vector<FabFactory<FAB> const*>& a_factory)
106{
107 define(a_geom, a_grids, a_dmap, a_overset_mask, a_info, a_factory);
108}
109
110template <typename MF>
111void
113 const Vector<BoxArray>& a_grids,
114 const Vector<DistributionMapping>& a_dmap,
115 const LPInfo& a_info,
116 const Vector<FabFactory<FAB> const*>& a_factory)
117{
118 BL_PROFILE("MLPoisson::define()");
119 MLCellABecLapT<MF>::define(a_geom, a_grids, a_dmap, a_info, a_factory);
120}
121
122template <typename MF>
123void
125 const Vector<BoxArray>& a_grids,
126 const Vector<DistributionMapping>& a_dmap,
127 const Vector<iMultiFab const*>& a_overset_mask,
128 const LPInfo& a_info,
129 const Vector<FabFactory<FAB> const*>& a_factory)
130{
131 BL_PROFILE("MLPoisson::define(overset)");
132 MLCellABecLapT<MF>::define(a_geom, a_grids, a_dmap, a_overset_mask, a_info, a_factory);
133}
134
135template <typename MF>
136MLPoissonT<MF>::~MLPoissonT () = default;
137
138template <typename MF>
139void
141{
142 BL_PROFILE("MLPoisson::prepareForSolve()");
143
145
146 m_is_singular.clear();
147 m_is_singular.resize(this->m_num_amr_levels, false);
148 auto itlo = std::find(this->m_lobc[0].begin(), this->m_lobc[0].end(), BCType::Dirichlet);
149 auto ithi = std::find(this->m_hibc[0].begin(), this->m_hibc[0].end(), BCType::Dirichlet);
150 if (itlo == this->m_lobc[0].end() && ithi == this->m_hibc[0].end())
151 { // No Dirichlet
152 for (int alev = 0; alev < this->m_num_amr_levels; ++alev)
153 {
154 // For now this assumes that overset regions are treated as Dirichlet bc's
155 if (this->m_domain_covered[alev] && !this->m_overset_mask[alev][0])
156 {
157 m_is_singular[alev] = true;
158 }
159 }
160 }
161 if (!m_is_singular[0] && this->m_needs_coarse_data_for_bc &&
162 this->m_coarse_fine_bc_type == LinOpBCType::Neumann)
163 {
164 AMREX_ASSERT(this->m_overset_mask[0][0] == nullptr);
165 auto bbox = this->m_grids[0][0].minimalBox();
166 for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
167 if (this->m_lobc[0][idim] == LinOpBCType::Dirichlet) {
168 bbox.growLo(idim,1);
169 }
170 if (this->m_hibc[0][idim] == LinOpBCType::Dirichlet) {
171 bbox.growHi(idim,1);
172 }
173 }
174 if (this->m_geom[0][0].Domain().contains(bbox)) {
175 m_is_singular[0] = true;
176 }
177 }
178}
179
180template <typename MF>
181void
182MLPoissonT<MF>::Fapply (int amrlev, int mglev, MF& out, const MF& in) const
183{
184 BL_PROFILE("MLPoisson::Fapply()");
185
186 const Real* dxinv = this->m_geom[amrlev][mglev].InvCellSize();
187
188 AMREX_D_TERM(const RT dhx = RT(dxinv[0]*dxinv[0]);,
189 const RT dhy = RT(dxinv[1]*dxinv[1]);,
190 const RT dhz = RT(dxinv[2]*dxinv[2]););
191
192#if (AMREX_SPACEDIM == 3)
193 RT dh0 = this->get_d0(dhx, dhy, dhz);
194 RT dh1 = this->get_d1(dhx, dhy, dhz);
195#endif
196
197#if (AMREX_SPACEDIM < 3)
198 const RT dx = RT(this->m_geom[amrlev][mglev].CellSize(0));
199 const RT probxlo = RT(this->m_geom[amrlev][mglev].ProbLo(0));
200#endif
201
202#ifdef AMREX_USE_GPU
203 if (Gpu::inLaunchRegion() && out.isFusingCandidate() && !this->hasHiddenDimension()) {
204 auto const& xma = in.const_arrays();
205 auto const& yma = out.arrays();
206 if (this->m_overset_mask[amrlev][mglev]) {
208 const auto& osmma = this->m_overset_mask[amrlev][mglev]->const_arrays();
209 ParallelFor(out,
210 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
211 {
213 mlpoisson_adotx_os(AMREX_D_DECL(i,j,k), yma[box_no], xma[box_no], osmma[box_no],
214 AMREX_D_DECL(dhx,dhy,dhz));
215 });
216 } else {
217#if (AMREX_SPACEDIM < 3)
218 if (this->m_has_metric_term) {
219 ParallelFor(out,
220 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
221 {
223 mlpoisson_adotx_m(AMREX_D_DECL(i,j,k), yma[box_no], xma[box_no],
224 AMREX_D_DECL(dhx,dhy,dhz), dx, probxlo);
225 });
226 } else
227#endif
228 {
229 ParallelFor(out,
230 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
231 {
233 mlpoisson_adotx(AMREX_D_DECL(i,j,k), yma[box_no], xma[box_no],
234 AMREX_D_DECL(dhx,dhy,dhz));
235 });
236 }
237 }
238 if (!Gpu::inNoSyncRegion()) {
240 }
241 } else
242#endif
243 {
244#ifdef AMREX_USE_OMP
245#pragma omp parallel if (Gpu::notInLaunchRegion())
246#endif
247 for (MFIter mfi(out, TilingIfNotGPU()); mfi.isValid(); ++mfi)
248 {
249 const Box& bx = mfi.tilebox();
250 const auto& xfab = in.array(mfi);
251 const auto& yfab = out.array(mfi);
252
253 if (this->m_overset_mask[amrlev][mglev]) {
255 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
257 {
259 mlpoisson_adotx_os(AMREX_D_DECL(i,j,k), yfab, xfab, osm,
260 AMREX_D_DECL(dhx,dhy,dhz));
261 });
262 } else {
263#if (AMREX_SPACEDIM == 3)
264 if (this->hasHiddenDimension()) {
265 Box const& bx2d = this->compactify(bx);
266 const auto& xfab2d = this->compactify(xfab);
267 const auto& yfab2d = this->compactify(yfab);
269 {
271 TwoD::mlpoisson_adotx(i, j, yfab2d, xfab2d, dh0, dh1);
272 });
273 } else {
275 {
276 mlpoisson_adotx(i, j, k, yfab, xfab, dhx, dhy, dhz);
277 });
278 }
279#elif (AMREX_SPACEDIM == 2)
280 if (this->m_has_metric_term) {
282 {
284 mlpoisson_adotx_m(i, j, yfab, xfab, dhx, dhy, dx, probxlo);
285 });
286 } else {
288 {
290 mlpoisson_adotx(i, j, yfab, xfab, dhx, dhy);
291 });
292 }
293#elif (AMREX_SPACEDIM == 1)
294 if (this->m_has_metric_term) {
296 {
298 mlpoisson_adotx_m(i, yfab, xfab, dhx, dx, probxlo);
299 });
300 } else {
302 {
304 mlpoisson_adotx(i, yfab, xfab, dhx);
305 });
306 }
307#endif
308 }
309 }
310 }
311}
312
313template <typename MF>
314void
315MLPoissonT<MF>::normalize (int amrlev, int mglev, MF& mf) const
316{
317 amrex::ignore_unused(amrlev,mglev,mf);
318#if (AMREX_SPACEDIM != 3)
319 BL_PROFILE("MLPoisson::normalize()");
320
321 if (!this->m_has_metric_term) { return; }
322
323 const Real* dxinv = this->m_geom[amrlev][mglev].InvCellSize();
324 AMREX_D_TERM(const RT dhx = RT(dxinv[0]*dxinv[0]);,
325 const RT dhy = RT(dxinv[1]*dxinv[1]);,
326 const RT dhz = RT(dxinv[2]*dxinv[2]););
327 const RT dx = RT(this->m_geom[amrlev][mglev].CellSize(0));
328 const RT probxlo = RT(this->m_geom[amrlev][mglev].ProbLo(0));
329
330#ifdef AMREX_USE_GPU
331 if (Gpu::inLaunchRegion() && mf.isFusingCandidate()) {
332 auto const& ma = mf.arrays();
333 ParallelFor(mf,
334 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
335 {
336 mlpoisson_normalize(i,j,k, ma[box_no], AMREX_D_DECL(dhx,dhy,dhz), dx, probxlo);
337 });
338 if (!Gpu::inNoSyncRegion()) {
340 }
341 } else
342#endif
343 {
344#ifdef AMREX_USE_OMP
345#pragma omp parallel if (Gpu::notInLaunchRegion())
346#endif
347 for (MFIter mfi(mf, TilingIfNotGPU()); mfi.isValid(); ++mfi)
348 {
349 const Box& bx = mfi.tilebox();
350 const auto& fab = mf.array(mfi);
351
352#if (AMREX_SPACEDIM == 2)
354 {
355 mlpoisson_normalize(i,j,k, fab, dhx, dhy, dx, probxlo);
356 });
357#else
359 {
360 mlpoisson_normalize(i,j,k, fab, dhx, dx, probxlo);
361 });
362#endif
363 }
364 }
365#endif
366}
367
368template <typename MF>
369void
370MLPoissonT<MF>::Fsmooth (int amrlev, int mglev, MF& sol, const MF& rhs, int redblack) const
371{
372 BL_PROFILE("MLPoisson::Fsmooth()");
373
374 MF Ax;
375 if (! this->m_use_gauss_seidel) { // jacobi
376 Ax.define(sol.boxArray(), sol.DistributionMap(), sol.nComp(), 0);
377 Fapply(amrlev, mglev, Ax, sol);
378 }
379
380 const auto& undrrelxr = this->m_undrrelxr[amrlev][mglev];
381 const auto& maskvals = this->m_maskvals [amrlev][mglev];
382
383 OrientationIter oitr;
384
385 const auto& f0 = undrrelxr[oitr()]; ++oitr;
386 const auto& f1 = undrrelxr[oitr()]; ++oitr;
387#if (AMREX_SPACEDIM > 1)
388 const auto& f2 = undrrelxr[oitr()]; ++oitr;
389 const auto& f3 = undrrelxr[oitr()]; ++oitr;
390#if (AMREX_SPACEDIM > 2)
391 const auto& f4 = undrrelxr[oitr()]; ++oitr;
392 const auto& f5 = undrrelxr[oitr()]; ++oitr;
393#endif
394#endif
395
396 const MultiMask& mm0 = maskvals[0];
397 const MultiMask& mm1 = maskvals[1];
398#if (AMREX_SPACEDIM > 1)
399 const MultiMask& mm2 = maskvals[2];
400 const MultiMask& mm3 = maskvals[3];
401#if (AMREX_SPACEDIM > 2)
402 const MultiMask& mm4 = maskvals[4];
403 const MultiMask& mm5 = maskvals[5];
404#endif
405#endif
406
407 const Real* dxinv = this->m_geom[amrlev][mglev].InvCellSize();
408 AMREX_D_TERM(const RT dhx = RT(dxinv[0]*dxinv[0]);,
409 const RT dhy = RT(dxinv[1]*dxinv[1]);,
410 const RT dhz = RT(dxinv[2]*dxinv[2]););
411
412#if (AMREX_SPACEDIM == 3)
413 RT dh0 = RT(this->get_d0(dhx, dhy, dhz));
414 RT dh1 = RT(this->get_d1(dhx, dhy, dhz));
415#endif
416
417#if (AMREX_SPACEDIM < 3)
418 const RT dx = RT(this->m_geom[amrlev][mglev].CellSize(0));
419 const RT probxlo = RT(this->m_geom[amrlev][mglev].ProbLo(0));
420#endif
421
422 MFItInfo mfi_info;
423 if (Gpu::notInLaunchRegion()) { mfi_info.EnableTiling().SetDynamic(true); }
424
425#ifdef AMREX_USE_GPU
426 if (Gpu::inLaunchRegion() && sol.isFusingCandidate()
427 && ! this->hasHiddenDimension())
428 {
429 const auto& m0ma = mm0.const_arrays();
430 const auto& m1ma = mm1.const_arrays();
431#if (AMREX_SPACEDIM > 1)
432 const auto& m2ma = mm2.const_arrays();
433 const auto& m3ma = mm3.const_arrays();
434#if (AMREX_SPACEDIM > 2)
435 const auto& m4ma = mm4.const_arrays();
436 const auto& m5ma = mm5.const_arrays();
437#endif
438#endif
439
440 const auto& solnma = sol.arrays();
441 const auto& rhsma = rhs.const_arrays();
442
443 AMREX_ALWAYS_ASSERT(rhs.nGrowVect() == 0);
444
445 const auto& f0ma = f0.const_arrays();
446 const auto& f1ma = f1.const_arrays();
447#if (AMREX_SPACEDIM > 1)
448 const auto& f2ma = f2.const_arrays();
449 const auto& f3ma = f3.const_arrays();
450#if (AMREX_SPACEDIM > 2)
451 const auto& f4ma = f4.const_arrays();
452 const auto& f5ma = f5.const_arrays();
453#endif
454#endif
455
456 if (this->m_overset_mask[amrlev][mglev]) {
458 const auto& osmma = this->m_overset_mask[amrlev][mglev]->const_arrays();
459 if (this->m_use_gauss_seidel) {
460 ParallelFor(sol,
461 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
462 {
463 Box vbx(rhsma[box_no]);
464 mlpoisson_gsrb_os(i, j, k, solnma[box_no], rhsma[box_no],
465 osmma[box_no], AMREX_D_DECL(dhx, dhy, dhz),
466 f0ma[box_no], m0ma[box_no],
467 f1ma[box_no], m1ma[box_no],
468#if (AMREX_SPACEDIM > 1)
469 f2ma[box_no], m2ma[box_no],
470 f3ma[box_no], m3ma[box_no],
471#if (AMREX_SPACEDIM > 2)
472 f4ma[box_no], m4ma[box_no],
473 f5ma[box_no], m5ma[box_no],
474#endif
475#endif
476 vbx, redblack);
477 });
478 } else {
479 const auto& axma = Ax.const_arrays();
480 ParallelFor(sol,
481 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
482 {
483 Box vbx(rhsma[box_no]);
484 mlpoisson_jacobi_os(i, j, k, solnma[box_no], rhsma[box_no],
485 axma[box_no], osmma[box_no],
486 AMREX_D_DECL(dhx, dhy, dhz),
487 f0ma[box_no], m0ma[box_no],
488 f1ma[box_no], m1ma[box_no],
489#if (AMREX_SPACEDIM > 1)
490 f2ma[box_no], m2ma[box_no],
491 f3ma[box_no], m3ma[box_no],
492#if (AMREX_SPACEDIM > 2)
493 f4ma[box_no], m4ma[box_no],
494 f5ma[box_no], m5ma[box_no],
495#endif
496#endif
497 vbx);
498 });
499 }
500 }
501#if (AMREX_SPACEDIM < 3)
502 else if (this->m_has_metric_term) {
503 if (this->m_use_gauss_seidel) {
504 ParallelFor(sol,
505 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
506 {
507 Box vbx(rhsma[box_no]);
508 mlpoisson_gsrb_m(i, j, k, solnma[box_no], rhsma[box_no],
509 AMREX_D_DECL(dhx, dhy, dhz),
510 f0ma[box_no], m0ma[box_no],
511 f1ma[box_no], m1ma[box_no],
512#if (AMREX_SPACEDIM > 1)
513 f2ma[box_no], m2ma[box_no],
514 f3ma[box_no], m3ma[box_no],
515#endif
516 vbx, redblack,
517 dx, probxlo);
518 });
519 } else {
520 const auto& axma = Ax.const_arrays();
521 ParallelFor(sol,
522 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
523 {
524 Box vbx(rhsma[box_no]);
525 mlpoisson_jacobi_m(i, j, k, solnma[box_no], rhsma[box_no],
526 axma[box_no], AMREX_D_DECL(dhx, dhy, dhz),
527 f0ma[box_no], m0ma[box_no],
528 f1ma[box_no], m1ma[box_no],
529#if (AMREX_SPACEDIM > 1)
530 f2ma[box_no], m2ma[box_no],
531 f3ma[box_no], m3ma[box_no],
532#endif
533 vbx, dx, probxlo);
534 });
535 }
536 }
537#endif
538 else {
539 if (this->m_use_gauss_seidel) {
540 ParallelFor(sol,
541 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
542 {
543 Box vbx(rhsma[box_no]);
544 mlpoisson_gsrb(i, j, k, solnma[box_no], rhsma[box_no],
545 AMREX_D_DECL(dhx, dhy, dhz),
546 f0ma[box_no], m0ma[box_no],
547 f1ma[box_no], m1ma[box_no],
548#if (AMREX_SPACEDIM > 1)
549 f2ma[box_no], m2ma[box_no],
550 f3ma[box_no], m3ma[box_no],
551#if (AMREX_SPACEDIM > 2)
552 f4ma[box_no], m4ma[box_no],
553 f5ma[box_no], m5ma[box_no],
554#endif
555#endif
556 vbx, redblack);
557 });
558 } else {
559 const auto& axma = Ax.const_arrays();
560 ParallelFor(sol,
561 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
562 {
563 Box vbx(rhsma[box_no]);
564 mlpoisson_jacobi(i, j, k, solnma[box_no], rhsma[box_no],
565 axma[box_no], AMREX_D_DECL(dhx, dhy, dhz),
566 f0ma[box_no], m0ma[box_no],
567 f1ma[box_no], m1ma[box_no],
568#if (AMREX_SPACEDIM > 1)
569 f2ma[box_no], m2ma[box_no],
570 f3ma[box_no], m3ma[box_no],
571#if (AMREX_SPACEDIM > 2)
572 f4ma[box_no], m4ma[box_no],
573 f5ma[box_no], m5ma[box_no],
574#endif
575#endif
576 vbx);
577 });
578 }
579 }
580 } else
581#endif
582 {
583#ifdef AMREX_USE_OMP
584#pragma omp parallel if (Gpu::notInLaunchRegion())
585#endif
586 for (MFIter mfi(sol,mfi_info); mfi.isValid(); ++mfi)
587 {
588 const auto& m0 = mm0.array(mfi);
589 const auto& m1 = mm1.array(mfi);
590#if (AMREX_SPACEDIM > 1)
591 const auto& m2 = mm2.array(mfi);
592 const auto& m3 = mm3.array(mfi);
593#if (AMREX_SPACEDIM > 2)
594 const auto& m4 = mm4.array(mfi);
595 const auto& m5 = mm5.array(mfi);
596#endif
597#endif
598
599 const Box& tbx = mfi.tilebox();
600 const Box& vbx = mfi.validbox();
601 const auto& solnfab = sol.array(mfi);
602 const auto& rhsfab = rhs.array(mfi);
603
604 const auto& f0fab = f0.array(mfi);
605 const auto& f1fab = f1.array(mfi);
606#if (AMREX_SPACEDIM > 1)
607 const auto& f2fab = f2.array(mfi);
608 const auto& f3fab = f3.array(mfi);
609#if (AMREX_SPACEDIM > 2)
610 const auto& f4fab = f4.array(mfi);
611 const auto& f5fab = f5.array(mfi);
612#endif
613#endif
614
615#if (AMREX_SPACEDIM == 1)
616 if (this->m_overset_mask[amrlev][mglev]) {
618 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
619 if (this->m_use_gauss_seidel) {
621 {
622 mlpoisson_gsrb_os(i, j, k, solnfab, rhsfab, osm, dhx,
623 f0fab, m0,
624 f1fab, m1,
625 vbx, redblack);
626 });
627 } else {
628 const auto& axfab = Ax.const_array(mfi);
630 {
631 mlpoisson_jacobi_os(i, j, k, solnfab, rhsfab, axfab,
632 osm, dhx,
633 f0fab, m0,
634 f1fab, m1,
635 vbx);
636 });
637 }
638 } else if (this->m_has_metric_term) {
639 if (this->m_use_gauss_seidel) {
641 {
642 mlpoisson_gsrb_m(i, j, k, solnfab, rhsfab, dhx,
643 f0fab, m0,
644 f1fab, m1,
645 vbx, redblack,
646 dx, probxlo);
647 });
648 } else {
649 const auto& axfab = Ax.const_array(mfi);
651 {
652 mlpoisson_jacobi_m(i, j, k, solnfab, rhsfab, axfab, dhx,
653 f0fab, m0,
654 f1fab, m1,
655 vbx, dx, probxlo);
656 });
657 }
658 } else {
659 if (this->m_use_gauss_seidel) {
661 {
662 mlpoisson_gsrb(i, j, k, solnfab, rhsfab, dhx,
663 f0fab, m0,
664 f1fab, m1,
665 vbx, redblack);
666 });
667 } else {
668 const auto& axfab = Ax.const_array(mfi);
670 {
671 mlpoisson_jacobi(i, j, k, solnfab, rhsfab, axfab, dhx,
672 f0fab, m0,
673 f1fab, m1,
674 vbx);
675 });
676 }
677 }
678#endif
679
680#if (AMREX_SPACEDIM == 2)
681 if (this->m_overset_mask[amrlev][mglev]) {
683 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
684 if (this->m_use_gauss_seidel) {
686 {
687 mlpoisson_gsrb_os(i, j, k, solnfab, rhsfab, osm, dhx, dhy,
688 f0fab, m0,
689 f1fab, m1,
690 f2fab, m2,
691 f3fab, m3,
692 vbx, redblack);
693 });
694 } else {
695 const auto& axfab = Ax.const_array(mfi);
697 {
698 mlpoisson_jacobi_os(i, j, k, solnfab, rhsfab, axfab,
699 osm, dhx, dhy,
700 f0fab, m0,
701 f1fab, m1,
702 f2fab, m2,
703 f3fab, m3,
704 vbx);
705 });
706 }
707 } else if (this->m_has_metric_term) {
708 if (this->m_use_gauss_seidel) {
710 {
711 mlpoisson_gsrb_m(i, j, k, solnfab, rhsfab, dhx, dhy,
712 f0fab, m0,
713 f1fab, m1,
714 f2fab, m2,
715 f3fab, m3,
716 vbx, redblack,
717 dx, probxlo);
718 });
719 } else {
720 const auto& axfab = Ax.const_array(mfi);
722 {
723 mlpoisson_jacobi_m(i, j, k, solnfab, rhsfab, axfab, dhx, dhy,
724 f0fab, m0,
725 f1fab, m1,
726 f2fab, m2,
727 f3fab, m3,
728 vbx, dx, probxlo);
729 });
730 }
731 } else {
732 if (this->m_use_gauss_seidel) {
734 {
735 mlpoisson_gsrb(i, j, k, solnfab, rhsfab, dhx, dhy,
736 f0fab, m0,
737 f1fab, m1,
738 f2fab, m2,
739 f3fab, m3,
740 vbx, redblack);
741 });
742 } else {
743 const auto& axfab = Ax.const_array(mfi);
745 {
746 mlpoisson_jacobi(i, j, k, solnfab, rhsfab, axfab, dhx, dhy,
747 f0fab, m0,
748 f1fab, m1,
749 f2fab, m2,
750 f3fab, m3,
751 vbx);
752 });
753 }
754 }
755#endif
756
757#if (AMREX_SPACEDIM == 3)
758 if (this->m_overset_mask[amrlev][mglev]) {
760 const auto& osm = this->m_overset_mask[amrlev][mglev]->const_array(mfi);
761 if (this->m_use_gauss_seidel) {
763 {
764 mlpoisson_gsrb_os(i, j, k, solnfab, rhsfab, osm, dhx, dhy, dhz,
765 f0fab, m0,
766 f1fab, m1,
767 f2fab, m2,
768 f3fab, m3,
769 f4fab, m4,
770 f5fab, m5,
771 vbx, redblack);
772 });
773 } else {
774 const auto& axfab = Ax.const_array(mfi);
776 {
777 mlpoisson_jacobi_os(i, j, k, solnfab, rhsfab, axfab,
778 osm, dhx, dhy, dhz,
779 f0fab, m0,
780 f1fab, m1,
781 f2fab, m2,
782 f3fab, m3,
783 f4fab, m4,
784 f5fab, m5,
785 vbx);
786 });
787 }
788 } else if (this->hasHiddenDimension()) {
789 Box const& tbx_2d = this->compactify(tbx);
790 Box const& vbx_2d = this->compactify(vbx);
791 const auto& solnfab_2d = this->compactify(solnfab);
792 const auto& rhsfab_2d = this->compactify(rhsfab);
793 const auto& f0fab_2d = this->compactify(this->get_d0(f0fab,f1fab,f2fab));
794 const auto& f1fab_2d = this->compactify(this->get_d1(f0fab,f1fab,f2fab));
795 const auto& f2fab_2d = this->compactify(this->get_d0(f3fab,f4fab,f5fab));
796 const auto& f3fab_2d = this->compactify(this->get_d1(f3fab,f4fab,f5fab));
797 const auto& m0_2d = this->compactify(this->get_d0(m0,m1,m2));
798 const auto& m1_2d = this->compactify(this->get_d1(m0,m1,m2));
799 const auto& m2_2d = this->compactify(this->get_d0(m3,m4,m5));
800 const auto& m3_2d = this->compactify(this->get_d1(m3,m4,m5));
801 if (this->m_use_gauss_seidel) {
802 AMREX_HOST_DEVICE_PARALLEL_FOR_3D ( tbx_2d, i, j, k,
803 {
804 TwoD::mlpoisson_gsrb(i, j, k, solnfab_2d, rhsfab_2d, dh0, dh1,
805 f0fab_2d, m0_2d,
806 f1fab_2d, m1_2d,
807 f2fab_2d, m2_2d,
808 f3fab_2d, m3_2d,
809 vbx_2d, redblack);
810 });
811 } else {
812 const auto& axfab = Ax.const_array(mfi);
813 const auto& axfab_2d = this->compactify(axfab);
814 AMREX_HOST_DEVICE_PARALLEL_FOR_3D ( tbx_2d, i, j, k,
815 {
816 TwoD::mlpoisson_jacobi(i, j, k, solnfab_2d, rhsfab_2d,
817 axfab_2d, dh0, dh1,
818 f0fab_2d, m0_2d,
819 f1fab_2d, m1_2d,
820 f2fab_2d, m2_2d,
821 f3fab_2d, m3_2d,
822 vbx_2d);
823 });
824 }
825 } else {
826 if (this->m_use_gauss_seidel) {
828 {
829 mlpoisson_gsrb(i, j, k, solnfab, rhsfab, dhx, dhy, dhz,
830 f0fab, m0,
831 f1fab, m1,
832 f2fab, m2,
833 f3fab, m3,
834 f4fab, m4,
835 f5fab, m5,
836 vbx, redblack);
837 });
838 } else {
839 const auto& axfab = Ax.const_array(mfi);
841 {
842 mlpoisson_jacobi(i, j, k, solnfab, rhsfab, axfab,
843 dhx, dhy, dhz,
844 f0fab, m0,
845 f1fab, m1,
846 f2fab, m2,
847 f3fab, m3,
848 f4fab, m4,
849 f5fab, m5,
850 vbx);
851 });
852 }
853 }
854#endif
855 }
856 }
857}
858
859template <typename MF>
860void
861MLPoissonT<MF>::FFlux (int amrlev, const MFIter& mfi,
862 const Array<FAB*,AMREX_SPACEDIM>& flux,
863 const FAB& sol, Location, const int face_only) const
864{
865 BL_PROFILE("MLPoisson::FFlux()");
866
867 const int mglev = 0;
868 const Box& box = mfi.tilebox();
869 const Real* dxinv = this->m_geom[amrlev][mglev].InvCellSize();
870
871 AMREX_D_TERM(const auto& fxarr = flux[0]->array();,
872 const auto& fyarr = flux[1]->array();,
873 const auto& fzarr = flux[2]->array(););
874 const auto& solarr = sol.array();
875
876#if (AMREX_SPACEDIM != 3)
877 const RT dx = RT(this->m_geom[amrlev][mglev].CellSize(0));
878 const RT probxlo = RT(this->m_geom[amrlev][mglev].ProbLo(0));
879#endif
880
881#if (AMREX_SPACEDIM == 3)
882 if (face_only) {
883 if (this->hiddenDirection() != 0) {
884 RT fac = RT(dxinv[0]);
885 Box blo = amrex::bdryLo(box, 0);
886 int blen = box.length(0);
888 {
889 mlpoisson_flux_xface(tbox, fxarr, solarr, fac, blen);
890 });
891 } else {
892 flux[0]->template setVal<RunOn::Device>(RT(0.0));
893 }
894 if (this->hiddenDirection() != 1) {
895 RT fac = RT(dxinv[1]);
896 Box blo = amrex::bdryLo(box, 1);
897 int blen = box.length(1);
899 {
900 mlpoisson_flux_yface(tbox, fyarr, solarr, fac, blen);
901 });
902 } else {
903 flux[1]->template setVal<RunOn::Device>(RT(0.0));
904 }
905 if (this->hiddenDirection() != 2) {
906 RT fac = RT(dxinv[2]);
907 Box blo = amrex::bdryLo(box, 2);
908 int blen = box.length(2);
910 {
911 mlpoisson_flux_zface(tbox, fzarr, solarr, fac, blen);
912 });
913 } else {
914 flux[2]->template setVal<RunOn::Device>(RT(0.0));
915 }
916 } else {
917 if (this->hiddenDirection() != 0) {
918 RT fac = RT(dxinv[0]);
919 Box bflux = amrex::surroundingNodes(box, 0);
921 {
922 mlpoisson_flux_x(tbox, fxarr, solarr, fac);
923 });
924 } else {
925 flux[0]->template setVal<RunOn::Device>(RT(0.0));
926 }
927 if (this->hiddenDirection() != 1) {
928 RT fac = RT(dxinv[1]);
929 Box bflux = amrex::surroundingNodes(box, 1);
931 {
932 mlpoisson_flux_y(tbox, fyarr, solarr, fac);
933 });
934 } else {
935 flux[1]->template setVal<RunOn::Device>(RT(0.0));
936 }
937 if (this->hiddenDirection() != 2) {
938 RT fac = RT(dxinv[2]);
939 Box bflux = amrex::surroundingNodes(box, 2);
941 {
942 mlpoisson_flux_z(tbox, fzarr, solarr, fac);
943 });
944 } else {
945 flux[2]->template setVal<RunOn::Device>(RT(0.0));
946 }
947 }
948#elif (AMREX_SPACEDIM == 2)
949 if (face_only) {
950 if (this->hiddenDirection() != 0) {
951 RT fac = RT(dxinv[0]);
952 Box blo = amrex::bdryLo(box, 0);
953 int blen = box.length(0);
954 if (this->m_has_metric_term) {
956 {
957 mlpoisson_flux_xface_m(tbox, fxarr, solarr, fac, blen, dx, probxlo);
958 });
959 } else {
961 {
962 mlpoisson_flux_xface(tbox, fxarr, solarr, fac, blen);
963 });
964 }
965 } else {
966 flux[0]->template setVal<RunOn::Device>(RT(0.0));
967 }
968 if (this->hiddenDirection() != 1) {
969 RT fac = RT(dxinv[1]);
970 Box blo = amrex::bdryLo(box, 1);
971 int blen = box.length(1);
972 if (this->m_has_metric_term) {
974 {
975 mlpoisson_flux_yface_m(tbox, fyarr, solarr, fac, blen, dx, probxlo);
976 });
977 } else {
979 {
980 mlpoisson_flux_yface(tbox, fyarr, solarr, fac, blen);
981 });
982 }
983 } else {
984 flux[1]->template setVal<RunOn::Device>(RT(0.0));
985 }
986 } else {
987 if (this->hiddenDirection() != 0) {
988 RT fac = RT(dxinv[0]);
989 Box bflux = amrex::surroundingNodes(box, 0);
990 if (this->m_has_metric_term) {
992 {
993 mlpoisson_flux_x_m(tbox, fxarr, solarr, fac, dx, probxlo);
994 });
995 } else {
997 {
998 mlpoisson_flux_x(tbox, fxarr, solarr, fac);
999 });
1000 }
1001 } else {
1002 flux[0]->template setVal<RunOn::Device>(RT(0.0));
1003 }
1004 if (this->hiddenDirection() != 1) {
1005 RT fac = RT(dxinv[1]);
1006 Box bflux = amrex::surroundingNodes(box, 1);
1007 if (this->m_has_metric_term) {
1009 {
1010 mlpoisson_flux_y_m(tbox, fyarr, solarr, fac, dx, probxlo);
1011 });
1012 } else {
1014 {
1015 mlpoisson_flux_y(tbox, fyarr, solarr, fac);
1016 });
1017 }
1018 } else {
1019 flux[1]->template setVal<RunOn::Device>(RT(0.0));
1020 }
1021 }
1022#else
1023 if (face_only) {
1024 RT fac = RT(dxinv[0]);
1025 Box blo = amrex::bdryLo(box, 0);
1026 int blen = box.length(0);
1027 if (this->m_has_metric_term) {
1029 {
1030 mlpoisson_flux_xface_m(tbox, fxarr, solarr, fac, blen, dx, probxlo);
1031 });
1032 } else {
1034 {
1035 mlpoisson_flux_xface(tbox, fxarr, solarr, fac, blen);
1036 });
1037 }
1038 } else {
1039 RT fac = RT(dxinv[0]);
1040 Box bflux = amrex::surroundingNodes(box, 0);
1041 if (this->m_has_metric_term) {
1043 {
1044 mlpoisson_flux_x_m(tbox, fxarr, solarr, fac, dx, probxlo);
1045 });
1046 } else {
1048 {
1049 mlpoisson_flux_x(tbox, fxarr, solarr, fac);
1050 });
1051 }
1052 }
1053#endif
1054}
1055
1056template <typename MF>
1057bool
1059{
1060 bool support = true;
1061 if (this->m_domain_covered[0]) { support = false; }
1062 if (this->doAgglomeration()) { support = false; }
1063 if (AMREX_SPACEDIM != 3) { support = false; }
1064 return support;
1065}
1066
1067template <typename MF>
1068std::unique_ptr<MLLinOpT<MF>>
1069MLPoissonT<MF>::makeNLinOp (int grid_size) const
1070{
1071 const Geometry& geom = this->m_geom[0].back();
1072 const BoxArray& ba = this->makeNGrids(grid_size);
1073
1075 {
1076 const std::vector<std::vector<int> >& sfc = DistributionMapping::makeSFC(ba);
1077 Vector<int> pmap(ba.size());
1079 const int nprocs = ParallelDescriptor::NProcs();
1080 for (int iproc = 0; iproc < nprocs; ++iproc) {
1081 for (int ibox : sfc[iproc]) {
1082 pmap[ibox] = iproc;
1083 }
1084 }
1085 dm.define(std::move(pmap));
1086 }
1087
1088 LPInfo minfo{};
1089 minfo.has_metric_term = this->info.has_metric_term;
1090
1091 std::unique_ptr<MLLinOpT<MF>> r{new MLALaplacianT<MF>({geom}, {ba}, {dm}, minfo)};
1092 auto nop = dynamic_cast<MLALaplacianT<MF>*>(r.get());
1093 if (!nop) {
1094 return nullptr;
1095 }
1096
1097 nop->m_parent = this;
1098
1099 nop->setMaxOrder(this->maxorder);
1100 nop->setVerbose(this->verbose);
1101
1102 nop->setDomainBC(this->m_lobc, this->m_hibc);
1103
1104 if (this->needsCoarseDataForBC())
1105 {
1106 const Real* dx0 = this->m_geom[0][0].CellSize();
1108 fac *= Real(0.5);
1109 RealVect cbloc {AMREX_D_DECL(dx0[0]*fac[0], dx0[1]*fac[1], dx0[2]*fac[2])};
1110 nop->setCoarseFineBCLocation(cbloc);
1111 }
1112
1113 nop->setScalars(1.0, -1.0);
1114
1115 const Real* dxinv = geom.InvCellSize();
1116 RT dxscale = RT(dxinv[0]);
1117#if (AMREX_SPACEDIM >= 2)
1118 dxscale = std::max(dxscale,RT(dxinv[1]));
1119#endif
1120#if (AMREX_SPACEDIM == 3)
1121 dxscale = std::max(dxscale,RT(dxinv[2]));
1122#endif
1123
1124 MF alpha(ba, dm, 1, 0);
1125 alpha.setVal(RT(1.e30)*dxscale*dxscale);
1126
1127 MF foo(this->m_grids[0].back(), this->m_dmap[0].back(), 1, 0, MFInfo().SetAlloc(false));
1128 const FabArrayBase::CPC& cpc = alpha.getCPC(IntVect(0),foo,IntVect(0),Periodicity::NonPeriodic());
1129 alpha.setVal(RT(0.0), cpc, 0, 1);
1130
1131 nop->setACoeffs(0, alpha);
1132
1133 return r;
1134}
1135
1136template <typename MF>
1137void
1138MLPoissonT<MF>::copyNSolveSolution (MF& dst, MF const& src) const
1139{
1140 dst.ParallelCopy(src);
1141}
1142
1143template <typename MF>
1144void
1146 MF const& phi)
1147{
1148 BL_PROFILE("MLPoisson::dpdn_faces()");
1149
1150 // We do not need to call applyBC because this function is used by the
1151 // OpenBC solver after solver has converged. That means the BC has been
1152 // filled to check the residual.
1153
1154 Box const& domain0 = this->m_geom[0][0].Domain();
1155 AMREX_D_TERM(const RT dxi = RT(this->m_geom[0][0].InvCellSize(0));,
1156 const RT dyi = RT(this->m_geom[0][0].InvCellSize(1));,
1157 const RT dzi = RT(this->m_geom[0][0].InvCellSize(2));)
1158
1159#ifdef AMREX_USE_OMP
1160#pragma omp parallel if (Gpu::notInLaunchRegion())
1161#endif
1162 for (MFIter mfi(phi); mfi.isValid(); ++mfi)
1163 {
1164 Box const& vbx = mfi.validbox();
1165 for (OrientationIter oit; oit.isValid(); ++oit) {
1166 Orientation face = oit();
1167 if (vbx[face] == domain0[face]) {
1168 int dir = face.coordDir();
1169 auto const& p = phi.const_array(mfi);
1170 auto const& gp = dpdn[dir]->array(mfi);
1171 Box const& b2d = amrex::bdryNode(vbx,face);
1172 if (dir == 0) {
1173 // because it's dphi/dn, not dphi/dx.
1174 RT fac = dxi * (face.isLow() ? RT(-1.0) : RT(1.));
1176 {
1177 gp(i,j,k) = fac * (p(i,j,k) - p(i-1,j,k));
1178 });
1179 }
1180#if (AMREX_SPACEDIM > 1)
1181 else if (dir == 1) {
1182 RT fac = dyi * (face.isLow() ? RT(-1.0) : RT(1.));
1184 {
1185 gp(i,j,k) = fac * (p(i,j,k) - p(i,j-1,k));
1186 });
1187 }
1188#if (AMREX_SPACEDIM > 2)
1189 else {
1190 RT fac = dzi * (face.isLow() ? RT(-1.0) : RT(1.));
1192 {
1193 gp(i,j,k) = fac * (p(i,j,k) - p(i,j,k-1));
1194 });
1195 }
1196#endif
1197#endif
1198 }
1199 }
1200 }
1201}
1202
1203extern template class MLPoissonT<MultiFab>;
1204
1206
1207}
1208
1209#endif
#define BL_PROFILE(a)
Definition AMReX_BLProfiler.H:551
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition AMReX_BLassert.H:50
#define AMREX_HOST_DEVICE_PARALLEL_FOR_3D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:110
#define AMREX_GPU_LAUNCH_HOST_DEVICE_LAMBDA_RANGE(TN, TI, block)
Definition AMReX_GpuLaunchMacrosC.nolint.H:4
#define AMREX_GPU_DEVICE
Definition AMReX_GpuQualifiers.H:18
#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
Long size() const noexcept
Return the number of boxes in the BoxArray.
Definition AMReX_BoxArray.H:597
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition AMReX_Box.H:146
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
void define(const BoxArray &boxes, int nprocs=ParallelDescriptor::NProcs())
Build mapping out of BoxArray over nprocs processors. You need to call this if you built your Distrib...
Definition AMReX_DistributionMapping.cpp:345
static DistributionMapping makeSFC(const MultiFab &weight, bool sort=true)
Definition AMReX_DistributionMapping.cpp:1764
Definition AMReX_FabFactory.H:50
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
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_MLALaplacian.H:14
Definition AMReX_MLCellABecLap.H:13
Vector< Vector< std::unique_ptr< iMultiFab > > > m_overset_mask
Definition AMReX_MLCellABecLap.H:85
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:249
BoxArray makeNGrids(int grid_size) const
Definition AMReX_MLCellLinOp.H:904
Vector< Vector< Array< MultiMask, 2 *AMREX_SPACEDIM > > > m_maskvals
Definition AMReX_MLCellLinOp.H:207
Vector< Vector< BndryRegisterT< MF > > > m_undrrelxr
Definition AMReX_MLCellLinOp.H:204
bool m_has_metric_term
Definition AMReX_MLCellLinOp.H:151
bool m_use_gauss_seidel
Definition AMReX_MLCellLinOp.H:213
T get_d0(T const &dx, T const &dy, T const &) const noexcept
Definition AMReX_MLLinOp.H:726
bool doAgglomeration() const noexcept
Definition AMReX_MLLinOp.H:681
Vector< Vector< BoxArray > > m_grids
Definition AMReX_MLLinOp.H:609
Vector< Vector< DistributionMapping > > m_dmap
Definition AMReX_MLLinOp.H:610
IntVect m_coarse_data_crse_ratio
Definition AMReX_MLLinOp.H:636
bool needsCoarseDataForBC() const noexcept
Needs coarse data for bc?
Definition AMReX_MLLinOp.H:177
bool hasHiddenDimension() const noexcept
Definition AMReX_MLLinOp.H:707
Vector< Array< BCType, AMREX_SPACEDIM > > m_lobc
Definition AMReX_MLLinOp.H:572
int hiddenDirection() const noexcept
Definition AMReX_MLLinOp.H:708
Vector< int > m_domain_covered
Definition AMReX_MLLinOp.H:612
const MLLinOpT< MF > * m_parent
Definition AMReX_MLLinOp.H:597
Vector< Vector< Geometry > > m_geom
first Vector is for amr level and second is mg level
Definition AMReX_MLLinOp.H:608
Box compactify(Box const &b) const noexcept
Definition AMReX_MLLinOp.H:1284
bool m_needs_coarse_data_for_bc
Definition AMReX_MLLinOp.H:634
int maxorder
Definition AMReX_MLLinOp.H:589
LPInfo info
Definition AMReX_MLLinOp.H:585
Vector< Array< BCType, AMREX_SPACEDIM > > m_hibc
Definition AMReX_MLLinOp.H:573
T get_d1(T const &, T const &dy, T const &dz) const noexcept
Definition AMReX_MLLinOp.H:736
LinOpBCType m_coarse_fine_bc_type
Definition AMReX_MLLinOp.H:635
int m_num_amr_levels
Definition AMReX_MLLinOp.H:593
Definition AMReX_MLPoisson.H:16
LinOpBCType BCType
Definition AMReX_MLPoisson.H:22
typename MF::value_type RT
Definition AMReX_MLPoisson.H:20
MLPoissonT< MF > & operator=(const MLPoissonT< MF > &)=delete
void copyNSolveSolution(MF &dst, MF const &src) const final
Definition AMReX_MLPoisson.H:1138
bool isBottomSingular() const final
Is the bottom of MG singular?
Definition AMReX_MLPoisson.H:59
void prepareForSolve() final
Definition AMReX_MLPoisson.H:140
void Fapply(int amrlev, int mglev, MF &out, const MF &in) const final
Definition AMReX_MLPoisson.H:182
void normalize(int amrlev, int mglev, MF &mf) const final
Divide mf by the diagonal component of the operator. Used by bicgstab.
Definition AMReX_MLPoisson.H:315
typename MF::fab_type FAB
Definition AMReX_MLPoisson.H:19
bool isSingular(int amrlev) const final
Is it singular on given AMR level?
Definition AMReX_MLPoisson.H:58
~MLPoissonT() override
MLPoissonT(const MLPoissonT< MF > &)=delete
bool supportNSolve() const final
Definition AMReX_MLPoisson.H:1058
MLPoissonT(MLPoissonT< MF > &&)=delete
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_MLPoisson.H:112
RT getAScalar() const final
Definition AMReX_MLPoisson.H:68
RT getBScalar() const final
Definition AMReX_MLPoisson.H:69
void get_dpdn_on_domain_faces(Array< MF *, AMREX_SPACEDIM > const &dpdn, MF const &phi)
Compute dphi/dn on domain faces after the solver has converged.
Definition AMReX_MLPoisson.H:1145
Array< MF const *, AMREX_SPACEDIM > getBCoeffs(int, int) const final
Definition AMReX_MLPoisson.H:71
MF const * getACoeffs(int, int) const final
Definition AMReX_MLPoisson.H:70
MLPoissonT()=default
Vector< int > m_is_singular
Definition AMReX_MLPoisson.H:86
void FFlux(int amrlev, const MFIter &mfi, const Array< FAB *, AMREX_SPACEDIM > &flux, const FAB &sol, Location loc, int face_only=0) const final
Definition AMReX_MLPoisson.H:861
typename MLLinOpT< MF >::Location Location
Definition AMReX_MLPoisson.H:23
void Fsmooth(int amrlev, int mglev, MF &sol, const MF &rhs, int redblack) const final
Definition AMReX_MLPoisson.H:370
std::unique_ptr< MLLinOpT< MF > > makeNLinOp(int grid_size) const final
Definition AMReX_MLPoisson.H:1069
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
AMREX_GPU_HOST_DEVICE bool isValid() const noexcept
Is the iterator valid?
Definition AMReX_Orientation.H:156
Encapsulation of the Orientation of the Faces of a Box.
Definition AMReX_Orientation.H:29
AMREX_GPU_HOST_DEVICE int coordDir() const noexcept
Returns the coordinate direction.
Definition AMReX_Orientation.H:83
AMREX_GPU_HOST_DEVICE bool isLow() const noexcept
Returns true if Orientation is low.
Definition AMReX_Orientation.H:89
static const Periodicity & NonPeriodic() noexcept
Definition AMReX_Periodicity.cpp:52
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
void streamSynchronize() noexcept
Definition AMReX_GpuDevice.H:237
bool inLaunchRegion() noexcept
Definition AMReX_GpuControl.H:86
bool notInLaunchRegion() noexcept
Definition AMReX_GpuControl.H:87
bool inNoSyncRegion() noexcept
Definition AMReX_GpuControl.H:146
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition AMReX_ParallelContext.H:70
MPI_Comm Communicator() noexcept
Definition AMReX_ParallelDescriptor.H:210
int NProcs() noexcept
return the number of MPI ranks local to the current Parallel Context
Definition AMReX_ParallelDescriptor.H:243
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi(int i, int j, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T dhx, T dhy, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Array4< T const > const &f2, Array4< int const > const &m2, Array4< T const > const &f3, Array4< int const > const &m3, Box const &vbox) noexcept
Definition AMReX_MLPoisson_2D_K.H:310
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx(int i, int j, Array4< T > const &y, Array4< T const > const &x, T dhx, T dhy) noexcept
Definition AMReX_MLPoisson_2D_K.H:13
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb(int i, int j, int, Array4< T > const &phi, Array4< T const > const &rhs, T dhx, T dhy, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Array4< T const > const &f2, Array4< int const > const &m2, Array4< T const > const &f3, Array4< int const > const &m3, Box const &vbox, int redblack) noexcept
Definition AMReX_MLPoisson_2D_K.H:197
Definition AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_yface(Box const &box, Array4< T > const &fy, Array4< T const > const &sol, T dyinv, int ylen) noexcept
Definition AMReX_MLPoisson_3D_K.H:90
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox) noexcept
Definition AMReX_MLPoisson_1D_K.H:193
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
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx_m(int i, Array4< T > const &y, Array4< T const > const &x, T dhx, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:32
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb_m(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:162
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 > bdryNode(const BoxND< dim > &b, Orientation face, int len=1) noexcept
Similar to bdryLo and bdryHi except that it operates on the given face of box b.
Definition AMReX_Box.H:1549
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_x(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv) noexcept
Definition AMReX_MLPoisson_1D_K.H:43
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx(int i, Array4< T > const &y, Array4< T const > const &x, T dhx) noexcept
Definition AMReX_MLPoisson_1D_K.H:9
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_y(Box const &box, Array4< T > const &fy, Array4< T const > const &sol, T dyinv) noexcept
Definition AMReX_MLPoisson_3D_K.H:72
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_xface_m(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, int xlen, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:86
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi_os(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, Array4< int const > const &osm, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox) noexcept
Definition AMReX_MLPoisson_1D_K.H:216
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_x_m(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:57
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack) noexcept
Definition AMReX_MLPoisson_1D_K.H:102
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 mlpoisson_flux_xface(Box const &box, Array4< T > const &fx, Array4< T const > const &sol, T dxinv, int xlen) noexcept
Definition AMReX_MLPoisson_1D_K.H:73
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 ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition AMReX.H:127
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_normalize(int i, int, int, Array4< T > const &x, T dhx, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:268
bool TilingIfNotGPU() noexcept
Definition AMReX_MFIter.H:12
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_gsrb_os(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< int const > const &osm, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, int redblack) noexcept
Definition AMReX_MLPoisson_1D_K.H:130
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition AMReX_Box.H:1881
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_zface(Box const &box, Array4< T > const &fz, Array4< T const > const &sol, T dzinv, int zlen) noexcept
Definition AMReX_MLPoisson_3D_K.H:130
int verbose
Definition AMReX_DistributionMapping.cpp:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_adotx_os(int i, Array4< T > const &y, Array4< T const > const &x, Array4< int const > const &osm, T dhx) noexcept
Definition AMReX_MLPoisson_1D_K.H:18
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_flux_z(Box const &box, Array4< T > const &fz, Array4< T const > const &sol, T dzinv) noexcept
Definition AMReX_MLPoisson_3D_K.H:112
std::array< T, N > Array
Definition AMReX_Array.H:24
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlpoisson_jacobi_m(int i, int, int, Array4< T > const &phi, Array4< T const > const &rhs, Array4< T const > const &Ax, T dhx, Array4< T const > const &f0, Array4< int const > const &m0, Array4< T const > const &f1, Array4< int const > const &m1, Box const &vbox, T dx, T probxlo) noexcept
Definition AMReX_MLPoisson_1D_K.H:242
parallel copy or add
Definition AMReX_FabArrayBase.H:536
Definition AMReX_MLLinOp.H:35
bool has_metric_term
Definition AMReX_MLLinOp.H:43
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