7 #include <AMReX_MultiFabUtil.H>
8 #include <AMReX_MLPoisson.H>
9 #include <AMReX_MLMG.H>
20 unsigned short baseLevel,
21 unsigned short finestLevel,
24 for (
int i = 0; i <= finestLevel; ++i) {
25 phi[i]->setVal(0.0, 1);
26 for (
int j = 0; j < AMREX_SPACEDIM; ++j) {
27 efield[i][j]->setVal(0.0, 1);
31 this->
mlmg_m(rho, phi, efield, baseLevel, finestLevel);
66 os <<
"* ************* A M R e X P o i s s o n S o l v e r ************************************ " <<
endl
67 <<
"* relative tolerance " <<
reltol_m <<
'\n'
69 <<
"* ******************************************************************** " <<
endl;
84 info.setAgglomeration(
true);
85 info.setConsolidation(
true);
86 info.setMaxCoarseningLevel(10);
90 int nlevels = finestLevel - baseLevel + 1;
97 for (
int ilev = 0; ilev < nlevels; ++ilev) {
98 geom_p[ilev] = geom[ilev + baseLevel];
99 ba[ilev] = rho[ilev + baseLevel]->boxArray();
100 dm[ilev] = rho[ilev + baseLevel]->DistributionMap();
101 rho_p[ilev] = rho[ilev + baseLevel].get();
102 phi_p[ilev] = phi[ilev + baseLevel].get();
105 amrex::MLPoisson mlpoisson(geom_p, ba, dm, info);
107 mlpoisson.setMaxOrder(3);
109 mlpoisson.setDomainBC({AMREX_D_DECL(amrex::LinOpBCType::Dirichlet,
110 amrex::LinOpBCType::Dirichlet,
111 amrex::LinOpBCType::Dirichlet)},
112 {AMREX_D_DECL(amrex::LinOpBCType::Dirichlet,
113 amrex::LinOpBCType::Dirichlet,
114 amrex::LinOpBCType::Dirichlet)});
116 for (
int ilev = 0; ilev < nlevels; ++ilev) {
117 mlpoisson.setLevelBC(ilev, phi[ilev].
get());
120 amrex::MLMG mlmg(mlpoisson);
121 mlmg.setMaxIter(200);
122 mlmg.setMaxFmgIter(0);
124 mlmg.setCGVerbose(0);
131 amrex::Vector<std::array<AmrField_t*, AMREX_SPACEDIM> > grad(nlevels);
133 for (
int i = 0; i < nlevels; ++i) {
134 for (
int j = 0; j < AMREX_SPACEDIM; ++j) {
135 grad[i][j] = efield[i][j].get();
139 mlmg.getGradSolution(grad);
141 for (
int i = 0; i < nlevels; ++i) {
142 for (
int j = 0; j < AMREX_SPACEDIM; ++j) {
143 efield[i][j]->mult(-1.0, 0, 1, 1);
double getZRangeMax(unsigned short level=0)
amrex::Vector< const AmrBoxLib::AmrField_t * > const_AmrFieldContainer_pt
constexpr double e
The value of .
AmrBoxLib::AmrGeomContainer_t GeomContainer_t
double getXRangeMax(unsigned short level=0)
double getYRangeMin(unsigned short level=0)
AmrBoxLib::AmrScalarFieldContainer_t AmrScalarFieldContainer_t
MLPoissonSolver(AmrBoxLib *itsAmrObject_p)
double getXRangeMin(unsigned short level=0)
amrex::Vector< AmrBoxLib::AmrField_t * > AmrFieldContainer_pt
double reltol_m
Relative tolearance for solver.
Inform & print(Inform &os) const
AmrBoxLib::AmrGridContainer_t AmrGridContainer_t
AmrBoxLib::AmrProcMapContainer_t AmrProcMapContainer_t
AmrBoxLib::AmrVectorFieldContainer_t AmrVectorFieldContainer_t
double getZRangeMin(unsigned short level=0)
AmrBoxLib * itsAmrObject_mp
double getYRangeMax(unsigned short level=0)
void solve(AmrScalarFieldContainer_t &rho, AmrScalarFieldContainer_t &phi, AmrVectorFieldContainer_t &efield, unsigned short baseLevel, unsigned short finestLevel, bool prevAsGuess=true)
double abstol_m
Absolute tolerance for solver.
Inform & endl(Inform &inf)
void mlmg_m(AmrScalarFieldContainer_t &rho, AmrScalarFieldContainer_t &phi, AmrVectorFieldContainer_t &efield, int baseLevel, int finestLevel)