Is there a conjugate gradient (CG) solver with MLMG as a preconditioner for multi-level MultiFab data?

[taiyaa]

I am currently working on some channel-like flow problem and noticed that if I use the default MLMG solver the norm of residual will grow gradually and the computation will blowup eventually. I solved the problem and created a demo by modifying a test problem in AMREX/Tests/LinearSolvers/CellEB/ and solved the Laplacian problem using gmres.solve() rather than mlmg.solve(). The gmres.solve() actually calls the mlmg as a preconditioner.

But my problem is the modification only works for single level MF data and it won't work for AMR mesh. Since it is quite common that people use Multigrid preconditioned CG/BiCGStab methods to solve flow problems, I would like to know if there are supports for such solvers?

The following is the code snippet I used in my problem:

MLEBABecLap mleb (geom, grids, dmap, info, amrex::GetVecOfConstPtrs(factory));
mleb.setMaxOrder(linop_maxorder);

mleb.setDomainBC(mlmg_lobc, mlmg_hibc);

// set phi value on cell centroid
// mleb.setPhiOnCentroid();

for (int ilev = 0; ilev <= max_level; ++ilev) {
    mleb.setLevelBC(ilev, &phi[ilev]);
}

mleb.setScalars(scalars[0], scalars[1]);

for (int ilev = 0; ilev <= max_level; ++ilev) {
    mleb.setACoeffs(ilev, acoef[ilev]);
    mleb.setBCoeffs(ilev, amrex::GetArrOfConstPtrs(bcoef[ilev]));
}

if (eb_is_dirichlet) {
    for (int ilev = 0; ilev <= max_level; ++ilev) {
        mleb.setEBDirichlet(ilev, phi_eb[ilev], bcoef_eb[ilev], phi_eb_type[ilev]);
    }
}

phi[0].setVal(10.0);


// solve use MLMG 
MLMG mlmg(mleb);
mlmg.setMaxIter(max_iter);
mlmg.setMaxFmgIter(max_fmg_iter);
mlmg.setBottomMaxIter(max_bottom_iter);
mlmg.setBottomTolerance(bottom_reltol);
mlmg.setVerbose(verbose);
mlmg.setBottomVerbose(bottom_verbose);
if (use_hypre) {
    mlmg.setBottomSolver(MLMG::BottomSolver::hypre);
} else if (use_petsc) {
    mlmg.setBottomSolver(MLMG::BottomSolver::petsc);
}
const Real tol_rel = reltol;
const Real tol_abs = 0.0;

// disable the mlmg solver
// mlmg.solve(amrex::GetVecOfPtrs(phi), amrex::GetVecOfConstPtrs(rhs), tol_rel, tol_abs);

// solve use GMRES
GMRESMLMGT<MultiFab> gmsolver(mlmg);
gmsolver.usePrecond(true);
gmsolver.setMaxIters(max_iter);
gmsolver.setVerbose(verbose);

// temporary array to check if the solution is converged
MultiFab phi_old;
phi_old.define(grids[0],dmap[0],2,0);

// assume that we have single amr level and solve the problem using GMRES solver
gmsolver.solve(phi[0], rhs[0], tol_rel, tol_abs);

[cgilet]

The default in IAMR is geometric multigrid with BiCGStab as the bottom solver. To change the bottom solver, for example to conjugate gradient, use mlmg.setBottomSolver(MLMG::BottomSolver::cg). Options for the bottom solver are in the documentation here: https://amrex-codes.github.io/amrex/docs_html/LinearSolvers.html?highlight=bottom+solver#parameters

Note that each linear solve is it's own thing, so the nodal projection is set separately from the MAC projection, which is separate from the diffusive/viscous solve. For the projections, you can change the bottom solver with the inputs file using, e.g., nodal_proj.bottom_solver = "cg" or mac_proj.bottom_solver = "bicgcg"

[taiyaa]

Sure. Let me clean up my code and will add the files in the attachment.

[taiyaa]

The demo is attached below.

CellEB_torus.zip

The result of phi is:

[WeiqunZhang]

I thought you had a test problem that mlmg failed to solve but gmres worked. So I was interested in testing with your problem. But the problem you provided does not appear to be the case. Both work just fine.

Note that in your test problem, the torus is very thin. So MLMG cannot coarsen the grids. It actually just uses bicgstab to solve the problem without any multigrid coarsening.

[WeiqunZhang]

Re: bottom solvers. Note that standard CG does not work with non-symmetric matrix. So you have to use bicg.

[taiyaa]

Thanks. I will pay more attention to the bottom solver for future cases. I might make some mistake when I configured the Torus case as I can't repeat the initial error message.
I do have a case that MLMG fails (but I am not able to test with GMRES as it involves two AMR levels). On the base grid the MLMG works but when I add one layer of AMR grid the computation will diverge. I reported the problem in this post(#174). I also attached my input file and the STL solid file if you are interested.

[WeiqunZhang]

Could you try AMReX-Codes/amrex#4271 to see if GMRES can help? If you notice any issues, let us know.

[taiyaa]

I tried the new GMRES solver (AMReX-Codes/amrex#4271). The problem in post (#174) is still there. It seems the preconditioned GMRES does not solve the convergence issue.