1 #define AMR_NO_SCALE false
4 template <
class MatrixType,
class VectorType>
7 const amrex::BoxArray& _grids,
8 const amrex::DistributionMapping& _dmap,
12 const Teuchos::RCP<comm_t>& comm,
13 const Teuchos::RCP<node_t>& node)
21 Bcrse_p(Teuchos::null),
22 Bfine_p(Teuchos::null),
26 residual_p(Teuchos::null),
27 error_p(Teuchos::null),
28 UnCovered_p(Teuchos::null),
33 for (
int j = 0; j < AMREX_SPACEDIM; ++j) {
34 G_p[j] = Teuchos::null;
36 nr_m[j] = _geom.Domain().length(j);
44 dx_m[j] = meshScaling[j] *
geom.CellSize(j);
61 template <
class MatrixType,
class VectorType>
64 map_p = Teuchos::null;
66 Anf_p = Teuchos::null;
69 Bcrse_p = Teuchos::null;
70 Bfine_p = Teuchos::null;
71 Awf_p = Teuchos::null;
73 for (
int j = 0; j < AMREX_SPACEDIM; ++j)
74 G_p[j] = Teuchos::null;
76 UnCovered_p = Teuchos::null;
78 rho_p = Teuchos::null;
79 phi_p = Teuchos::null;
80 residual_p = Teuchos::null;
81 error_p = Teuchos::null;
85 template <
class MatrixType,
class VectorType>
88 #if AMREX_SPACEDIM == 3
89 return iv[0] + (iv[1] + nr_m[1] * iv[2]) * nr_m[0];
91 return iv[0] + iv[1] * nr_m[0];
96 template <
class MatrixType,
class VectorType>
99 return bc_mp[0]->isBoundary(iv, &nr_m[0]);
103 template <
class MatrixType,
class VectorType>
108 bool applied =
false;
109 for (
int d = 0; d < AMREX_SPACEDIM; ++d) {
110 if ( bc_mp[d]->isBoundary(iv, d, &nr_m[0]) ) {
112 bc_mp[d]->apply(iv, d, map, value,
this, &nr_m[0]);
119 template <
class MatrixType,
class VectorType>
125 if ( fab(iv) != Mask::PHYSBNDRY )
128 bool applied =
false;
129 for (
int d = 0; d < AMREX_SPACEDIM; ++d) {
130 if ( bc_mp[d]->isBoundary(iv, d, &nr_m[0]) ) {
132 bc_mp[d]->apply(iv, d, map, value,
this, &nr_m[0]);
139 template <
class MatrixType,
class VectorType>
145 bc_mp[dir]->apply(iv, dir, map, value,
this, &nr_m[0]);
149 template <
class MatrixType,
class VectorType>
151 amrex::Periodicity period(
AmrIntVect_t(D_DECL(0, 0, 0)));
152 mask.reset(
new mask_t(grids, dmap, 1, 1));
153 mask->BuildMask(geom.Domain(), period,
154 Mask::COVERED, Mask::BNDRY,
155 Mask::PHYSBNDRY, Mask::INTERIOR);
156 mask->FillBoundary(period);
160 template <
class MatrixType,
class VectorType>
166 template <
class MatrixType,
class VectorType>
172 template <
class MatrixType,
class VectorType>
178 template <
class MatrixType,
class VectorType>
184 template <
class MatrixType,
class VectorType>
190 template <
class MatrixType,
class VectorType>
193 iv.allLT(
AmrIntVect_t(D_DECL(nr_m[0], nr_m[1], nr_m[2]))) );
197 template <
class MatrixType,
class VectorType>
199 const Teuchos::RCP<node_t>& node)
202 go_t localNumElements = 0;
204 Teuchos::Array<go_t> globalindices;
206 for (amrex::MFIter mfi(grids, dmap,
true); mfi.isValid(); ++mfi) {
207 const amrex::Box& tbx = mfi.tilebox();
208 const int* lo = tbx.loVect();
209 const int* hi = tbx.hiVect();
211 for (
int i = lo[0]; i <= hi[0]; ++i) {
212 for (
int j = lo[1]; j <= hi[1]; ++j) {
213 #if AMREX_SPACEDIM == 3
214 for (
int k = lo[2]; k <= hi[2]; ++k) {
220 globalindices.push_back(globalidx);
223 #if AMREX_SPACEDIM == 3
235 amrex::Box bx = grids.minimalBox();
236 const int* lo = bx.loVect();
243 go_t N = grids.numPts();
245 map_p = Teuchos::rcp(
new dmap_t(N, globalindices, baseIndex, comm, node) );
bool applyBoundary(const AmrIntVect_t &iv, umap_t &map, const scalar_t &value)
std::unordered_map< go_t, scalar_t > umap_t
Teuchos::RCP< matrix_t > G_p[AMREX_SPACEDIM]
gradient matrices in x, y, and z to compute electric field
amrex::IntVect AmrIntVect_t
void serialize(Param_t params, std::ostringstream &os)
serializes params using Boost archive
scalar_t invdx_m[AMREX_SPACEDIM]
inverse cell size in particle rest frame
const scalar_t * cellSize() const
amr::local_ordinal_t lo_t
bool isValid(const AmrIntVect_t &iv) const
boundary_t bc_mp[AMREX_SPACEDIM]
boundary conditions
amrex::FabArray< basefab_t > mask_t
amr::AmrIntVect_t AmrIntVect_t
void buildMap(const Teuchos::RCP< comm_t > &comm, const Teuchos::RCP< node_t > &node)
scalar_t dx_m[AMREX_SPACEDIM]
cell size in particle rest frame
std::shared_ptr< AmrBoundary< AmrMultiGridLevel< MatrixType, VectorType > > > boundary_t
Tpetra::Map< local_ordinal_t, global_ordinal_t, node_t > dmap_t
bool isBoundary(const AmrIntVect_t &iv) const
Teuchos::RCP< vector_t > residual_p
residual over all cells
amr::global_ordinal_t go_t
const AmrGeometry_t & geom
geometry of this problem
Teuchos::RCP< vector_t > error_p
error over all cells
amrex::BaseFab< int > basefab_t
go_t serialize(const AmrIntVect_t &iv) const
const AmrIntVect_t & refinement() const
Teuchos::RCP< dmap_t > map_p
Tpetra core map.
go_t nr_m[AMREX_SPACEDIM]
number of grid points
const scalar_t * invCellSize() const
amr::AmrGeometry_t AmrGeometry_t