Commit 3ee74cf1 authored by frey_m's avatar frey_m Committed by frey_m

Resolve "compiler errors with AMR / SAAMG enabled"

parent a1d2c53b
......@@ -51,6 +51,11 @@ add_compile_options (-Werror)
add_compile_options (-funroll-loops)
add_compile_options (-fstrict-aliasing)
# An error in Trilinos/12.18.1 causes
# -Waggressive-loop-optimizations with -O3 optimization.
# See issue 587 and corresponding merge request 415.
add_compile_options (-fno-aggressive-loop-optimizations)
option (USE_STATIC_LIBRARIES "Link with static libraries if available" ON)
if (NOT ${USE_STATIC_LIBRARIES})
add_compile_options (-fPIE -fPIC)
......@@ -106,6 +111,9 @@ if (CMAKE_Fortran_COMPILER)
if (${CMAKE_CXX_COMPILER_ID} STREQUAL "GNU")
set (Fortran_LIBS "-lgfortran")
endif ()
# CMake issue. Missing library libmpi_mpifh in MPI_Fortran_LINRARIES.
# See issue 587 and corresponding merge request 415.
set (Fortran_LIBS "${Fortran_LIBS} -lmpi_mpifh")
else ()
message (STATUS "No Fortran compiler found (optional)")
if (ENABLE_AMR)
......
......@@ -41,7 +41,7 @@ void Amesos2BottomSolver<Level>::solve(const Teuchos::RCP<mv_t>& x,
template <class Level>
void Amesos2BottomSolver<Level>::setOperator(const Teuchos::RCP<matrix_t>& A,
Level* level_p)
Level* /*level_p*/)
{
try {
solver_mp = Amesos2::create<matrix_t, mv_t>(solvertype_m, A);
......
......@@ -43,7 +43,7 @@ public:
umap_t& map,
const scalar_t& value,
Level* mglevel,
const go_t* nr);
const go_t* /*nr*/);
};
......@@ -53,7 +53,7 @@ void AmrDirichletBoundary<Level>::apply(const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& value,
Level* mglevel,
const go_t* nr)
const go_t* /*nr*/)
{
// find interior neighbour cell
AmrIntVect_t niv = iv;
......
......@@ -51,6 +51,7 @@ public:
/*!
* Get the stencil to interpolate a value from coarse to fine level
* @param iv is the fine cell where we want to have the interpolated value
* @param fab cells
* @param map with global matrix indices of coarse level cells and
* matrix entries of coarse level cells (coefficients)
* @param scale to apply to matrix values
......@@ -74,7 +75,7 @@ public:
* @param shift is either -1 or 1. If the refined coarse cell is on the left / lower / front
* side, shift is equal to -1, otherwise the interface is on the right / upper / back side
* and the value is 1.
* @param ba contains all coarse cells that got refined
* @param rfab contains all coarse cells that got refined
* @param riv is the fine cell at the interface
* @param mglevel used to get the global indices and refinement ratio among levels,
* and boundary values at physical domain, e.g. Dirichlet, open BC
......@@ -98,15 +99,15 @@ public:
* @param shift is either -1 or 1. If the refined coarse cell is on the left / lower / front
* side, shift is equal to -1, otherwise the interface is on the right / upper / back side
* and the value is 1.
* @param ba contains all coarse cells that got refined
* @param fab contains all coarse cells that got refined
* @param mglevel used to get the global indices and refinement ratio among levels,
* and boundary avlues at physical domain, e.g. Dirichlet, open BC
*/
virtual void fine(const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& fab,
Level* mglevel)
virtual void fine(const AmrIntVect_t& /*iv*/,
umap_t& /*map*/,
const scalar_t& /*scale*/,
lo_t /*dir*/, lo_t /*shift*/, const basefab_t& /*fab*/,
Level* /*mglevel*/)
{ };
/*!
......@@ -124,11 +125,11 @@ public:
* @param mglevel used to get the global indices and refinement ratio among levels,
* and boundary avlues at physical domain, e.g. Dirichlet, open BC
*/
virtual void fine(const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift,
Level* mglevel)
virtual void fine(const AmrIntVect_t& /*iv*/,
umap_t& /*map*/,
const scalar_t& /*scale*/,
lo_t /*dir*/, lo_t /*shift*/,
Level* /*mglevel*/)
{ };
protected:
......
......@@ -93,11 +93,11 @@ AmrLagrangeInterpolater<Level>::AmrLagrangeInterpolater(Order order)
template <class Level>
void AmrLagrangeInterpolater<Level>::stencil(
const AmrIntVect_t& iv,
const basefab_t& fab,
typename Level::umap_t& map,
const typename Level::scalar_t& scale,
Level* mglevel)
const AmrIntVect_t& /*iv*/,
const basefab_t& /*fab*/,
typename Level::umap_t& /*map*/,
const typename Level::scalar_t& /*scale*/,
Level* /*mglevel*/)
{
}
......@@ -198,7 +198,7 @@ void AmrLagrangeInterpolater<Level>::crseLinear_m(
const AmrIntVect_t& iv,
typename Level::umap_t& map,
const typename Level::scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& rfab,
lo_t dir, lo_t /*shift*/, const basefab_t& rfab,
const AmrIntVect_t& riv,
Level* mglevel)
{
......
......@@ -30,6 +30,7 @@
#include "AbstractObjects/OpalData.h"
#include "Utilities/OpalException.h"
#include "Utilities/Timer.h"
#include "Utilities/Util.h"
#include <AMReX_ParallelDescriptor.H>
#include <boost/filesystem.hpp>
......@@ -67,8 +68,6 @@ AmrMultiGrid::AmrMultiGrid(AmrBoxLib* itsAmrObject_p,
, fname_m(OpalData::getInstance()->getInputBasename() + std::string(".solver"))
, flag_m(std::ios::out)
{
node_mp = KokkosClassic::Details::getNode<amr::node_t>(); //KokkosClassic::DefaultNode::getDefaultNode();
#if AMR_MG_TIMER
this->initTimer_m();
#endif
......@@ -93,7 +92,7 @@ AmrMultiGrid::AmrMultiGrid(AmrBoxLib* itsAmrObject_p,
// preconditioner
const Preconditioner precond = this->convertToEnumPreconditioner_m(prec);
this->initPrec_m(precond, rebalance, reuse);
this->initPrec_m(precond, reuse);
// base level solver
const BaseSolver solver = this->convertToEnumBaseSolver_m(bsolver);
......@@ -164,10 +163,6 @@ void AmrMultiGrid::solve(AmrScalarFieldContainer_t &rho,
void AmrMultiGrid::setNumberOfSweeps(const std::size_t& nSweeps) {
if ( nSweeps < 0 )
throw OpalException("AmrMultiGrid::setNumberOfSweeps()",
"The number of smoothing sweeps needs to be non-negative!");
nSweeps_m = nSweeps;
}
......@@ -258,8 +253,7 @@ void AmrMultiGrid::initLevels_m(const amrex::Vector<AmrField_u>& rho,
geom[ilev],
rr,
bc_m,
comm_mp,
node_mp));
comm_mp));
} else {
mglevel_m[lev]->buildLevelMask();
}
......@@ -844,10 +838,10 @@ void AmrMultiGrid::buildMultiLevel_m(const amrex::Vector<AmrField_u>& rho,
cfab, invdx2);
this->buildFineBoundaryMatrix_m(lev, gidx, iv,
mfab, rfab, cfab);
mfab, rfab);
this->buildCompositePoissonMatrix_m(lev, gidx, iv, mfab,
rfab, cfab, invdx2);
rfab, invdx2);
if (lev > lbase_m || (lev == lbase_m && !solver_mp->hasOperator())) {
this->buildNoFinePoissonMatrix_m(lev, gidx, iv, mfab, invdx2);
......@@ -1119,7 +1113,6 @@ void AmrMultiGrid::buildCompositePoissonMatrix_m(const lo_t& level,
const AmrIntVect_t& iv,
const basefab_t& mfab,
const basefab_t& rfab,
const basefab_t& cfab,
const scalar_t* invdx2)
{
/*
......@@ -1431,8 +1424,7 @@ void AmrMultiGrid::buildFineBoundaryMatrix_m(const lo_t& level,
const go_t& gidx,
const AmrIntVect_t& iv,
const basefab_t& mfab,
const basefab_t& rfab,
const basefab_t& cfab)
const basefab_t& rfab)
{
/* fine: level + 1
* coarse (this): level
......@@ -1455,8 +1447,7 @@ void AmrMultiGrid::buildFineBoundaryMatrix_m(const lo_t& level,
auto fill = [&](umap_t& map,
D_DECL(int ii, int jj, int kk),
int* begin, int* end, int d,
const AmrIntVect_t& iv, int shift,
int sign)
const AmrIntVect_t& iv, int shift)
{
for (int iref = ii - begin[0]; iref <= ii + end[0]; ++iref) {
for (int jref = jj - begin[1]; jref <= jj + end[1]; ++jref) {
......@@ -1539,7 +1530,7 @@ void AmrMultiGrid::buildFineBoundaryMatrix_m(const lo_t& level,
// iterate over all fine cells at the interface
// start with lower cells --> cover coarse neighbour
// cell
fill(map, D_DECL(ii, jj, kk), &begin[0], &end[0], d, iv, shift, 1.0);
fill(map, D_DECL(ii, jj, kk), &begin[0], &end[0], d, iv, shift);
break;
}
case 1:
......@@ -1551,7 +1542,7 @@ void AmrMultiGrid::buildFineBoundaryMatrix_m(const lo_t& level,
#if AMREX_SPACEDIM == 3
int kk = covered[2] << 1; // refinemet in z
#endif
fill(map, D_DECL(ii, jj, kk), &begin[0], &end[0], d, iv, shift, 1.0);
fill(map, D_DECL(ii, jj, kk), &begin[0], &end[0], d, iv, shift);
break;
}
}
......@@ -1768,7 +1759,7 @@ void AmrMultiGrid::smooth_m(const lo_t& level,
#endif
// base level has no smoother --> l - 1
smoother_m[level-1]->smooth(e, mglevel_m[level]->Anf_p, r);
smoother_m[level-1]->smooth(e, r);
#if AMR_MG_TIMER
IpplTimings::stopTimer(smoothTimer_m);
......@@ -1985,7 +1976,6 @@ void AmrMultiGrid::initBaseSolver_m(const BaseSolver& solver,
void AmrMultiGrid::initPrec_m(const Preconditioner& prec,
const bool& rebalance,
const std::string& reuse)
{
switch ( prec ) {
......@@ -2001,7 +1991,7 @@ void AmrMultiGrid::initPrec_m(const Preconditioner& prec,
case Preconditioner::SA:
{
std::string muelu = MueLuPreconditioner_t::convertToMueLuReuseOption(reuse);
prec_mp.reset( new MueLuPreconditioner_t(rebalance, muelu) );
prec_mp.reset( new MueLuPreconditioner_t(muelu) );
break;
}
case Preconditioner::NONE:
......
......@@ -410,7 +410,6 @@ private:
const AmrIntVect_t& iv,
const basefab_t& mfab,
const basefab_t& rfab,
const basefab_t& cfab,
const scalar_t* invdx2);
/*!
......@@ -473,18 +472,17 @@ private:
* x^{(l)} = B_{fine}\cdot x^{(l+1)}
* \f]
* Dirichlet boundary condition. Flux matching.
* @param iv is the current cell
* @param cells all coarse cells that are at the crse-fine interface but are
* not refined
* @param crse_fine_ba coarse cells that got refined
* @param level the finest level is omitted
* @param gidx the global index
* @param iv is the current cell
* @param mfab is the mask (internal cell, boundary cell, ...) of that level
* @param rfab is the mask between levels
*/
void buildFineBoundaryMatrix_m(const lo_t& level,
const go_t& gidx,
const AmrIntVect_t& iv,
const basefab_t& mfab,
const basefab_t& rfab,
const basefab_t& cfab);
const basefab_t& rfab);
/*!
* Copy data from AMReX to Trilinos
......@@ -601,11 +599,9 @@ private:
/*!
* Instantiate a preconditioner for the bottom solver
* @param precond type
* @param rebalance preconditioner (SA only)
* @param reuse types of SA hierarchy
*/
void initPrec_m(const Preconditioner& prec,
const bool& rebalance,
const std::string& reuse);
/*!
......@@ -665,7 +661,6 @@ private:
private:
Teuchos::RCP<comm_t> comm_mp; ///< communicator
Teuchos::RCP<amr::node_t> node_mp; ///< kokkos node
/// interpolater without coarse-fine interface
std::unique_ptr<AmrInterpolater<AmrMultiGridLevel_t> > interp_mp;
......
......@@ -53,7 +53,6 @@ public:
typedef amr::comm_t comm_t;
typedef amr::dmap_t dmap_t;
typedef amr::node_t node_t;
typedef amr::global_ordinal_t go_t;
typedef amr::scalar_t scalar_t;
typedef amr::local_ordinal_t lo_t;
......@@ -96,7 +95,6 @@ public:
* @param rr refinement ratio
* @param bc physical boundaries (x, y, z)
* @param comm MPI communicator
* @param node Kokkos node type (Serial, OpenMP, CUDA)
*/
AmrMultiGridLevel(const Vector_t& meshScaling,
const amrex::BoxArray& _grids,
......@@ -104,8 +102,7 @@ public:
const AmrGeometry_t& _geom,
const AmrIntVect_t& rr,
const boundary_t* bc,
const Teuchos::RCP<comm_t>& comm,
const Teuchos::RCP<node_t>& node);
const Teuchos::RCP<comm_t>& comm);
~AmrMultiGridLevel();
......@@ -198,10 +195,8 @@ private:
/*!
* Build Tpetra::Map of this level
* @param comm MPI communicator
* @param node Kokkos node type
*/
void buildMap(const Teuchos::RCP<comm_t>& comm,
const Teuchos::RCP<node_t>& node);
void buildMap(const Teuchos::RCP<comm_t>& comm);
public:
const amrex::BoxArray& grids; ///< boxes of this level
......
......@@ -30,8 +30,7 @@ AmrMultiGridLevel<MatrixType,
const AmrGeometry_t& _geom,
const AmrIntVect_t& rr,
const boundary_t* bc,
const Teuchos::RCP<comm_t>& comm,
const Teuchos::RCP<node_t>& node)
const Teuchos::RCP<comm_t>& comm)
: grids(_grids),
dmap(_dmap),
geom(_geom),
......@@ -71,7 +70,7 @@ AmrMultiGridLevel<MatrixType,
this->buildLevelMask();
this->buildMap(comm, node);
this->buildMap(comm);
residual_p = Teuchos::rcp( new vector_t(map_p, false) );
......@@ -216,8 +215,7 @@ bool AmrMultiGridLevel<MatrixType, VectorType>::isValid(const AmrIntVect_t& iv)
template <class MatrixType, class VectorType>
void AmrMultiGridLevel<MatrixType, VectorType>::buildMap(const Teuchos::RCP<comm_t>& comm,
const Teuchos::RCP<node_t>& node)
void AmrMultiGridLevel<MatrixType, VectorType>::buildMap(const Teuchos::RCP<comm_t>& comm)
{
go_t localNumElements = 0;
......@@ -263,5 +261,5 @@ void AmrMultiGridLevel<MatrixType, VectorType>::buildMap(const Teuchos::RCP<comm
// numGlobalElements == N
go_t N = grids.numPts();
map_p = Teuchos::rcp( new dmap_t(N, globalindices, baseIndex, comm, node) );
map_p = Teuchos::rcp( new dmap_t(N, globalindices, baseIndex, comm) );
}
......@@ -42,12 +42,12 @@ void AmrPCInterpolater<Level>::stencil(
template <class Level>
void AmrPCInterpolater<Level>::coarse(
const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& rfab,
const AmrIntVect_t& riv,
Level* mglevel)
const AmrIntVect_t& /*iv*/,
umap_t& /*map*/,
const scalar_t& /*scale*/,
lo_t /*dir*/, lo_t /*shift*/, const basefab_t& /*rfab*/,
const AmrIntVect_t& /*riv*/,
Level* /*mglevel*/)
{
// do nothing
}
......@@ -58,7 +58,7 @@ void AmrPCInterpolater<Level>::fine(
const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& fab,
lo_t /*dir*/, lo_t /*shift*/, const basefab_t& fab,
Level* mglevel)
{
/*
......
......@@ -52,7 +52,6 @@ AmrSmoother::~AmrSmoother() {
void AmrSmoother::smooth(const Teuchos::RCP<vector_t>& x,
const Teuchos::RCP<matrix_t>& A,
const Teuchos::RCP<vector_t>& b)
{
prec_mp->apply(*b, *x, Teuchos::NO_TRANS,
......
......@@ -65,11 +65,9 @@ public:
/*!
* Perform one smoothing step
* @param x right-hand side
* @param A system matrix
* @param b right-hand side
*/
void smooth(const Teuchos::RCP<vector_t>& x,
const Teuchos::RCP<matrix_t>& A,
const Teuchos::RCP<vector_t>& b);
/*!
......
......@@ -116,12 +116,12 @@ void AmrTrilinearInterpolater<Level>::stencil(
template <class Level>
void AmrTrilinearInterpolater<Level>::coarse(
const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& rfab,
const AmrIntVect_t& riv,
Level* mglevel)
const AmrIntVect_t& /*iv*/,
umap_t& /*map*/,
const scalar_t& /*scale*/,
lo_t /*dir*/, lo_t /*shift*/, const basefab_t& /*rfab*/,
const AmrIntVect_t& /*riv*/,
Level* /*mglevel*/)
{
// do nothing
}
......@@ -132,7 +132,7 @@ void AmrTrilinearInterpolater<Level>::fine(
const AmrIntVect_t& iv,
umap_t& map,
const scalar_t& scale,
lo_t dir, lo_t shift, const basefab_t& fab,
lo_t /*dir*/, lo_t /*shift*/, const basefab_t& fab,
Level* mglevel)
{
/*
......
......@@ -72,7 +72,7 @@ void BelosBottomSolver<Level>::setOperator(const Teuchos::RCP<matrix_t>& A,
Level* level_p)
{
// make positive definite --> rhs has to change sign as well
A_mp = A->clone(A->getNode());
A_mp = Teuchos::rcp(new matrix_t(*A, Teuchos::Copy));
A_mp->resumeFill();
A_mp->scale(-1.0);
A_mp->fillComplete();
......
......@@ -30,7 +30,7 @@ Ifpack2Preconditioner<Level>::Ifpack2Preconditioner(Preconditioner prec)
template <class Level>
void Ifpack2Preconditioner<Level>::create(const Teuchos::RCP<amr::matrix_t>& A,
Level* level_p)
Level* /*level_p*/)
{
Ifpack2::Factory factory;
......
......@@ -52,8 +52,7 @@ public:
public:
MueLuPreconditioner(const bool& rebalance,
const std::string& reuse);
explicit MueLuPreconditioner(const std::string& reuse);
void create(const Teuchos::RCP<amr::matrix_t>& A, Level* level_p = nullptr);
......@@ -70,10 +69,6 @@ private:
Teuchos::ParameterList params_m;
Teuchos::RCP<precond_t> prec_mp;
Teuchos::RCP<amr::multivector_t> coords_mp;
const bool rebalance_m;
};
#include "MueLuPreconditioner.hpp"
......
......@@ -24,11 +24,8 @@
#include <MueLu_CreateTpetraPreconditioner.hpp>
template <class Level>
MueLuPreconditioner<Level>::MueLuPreconditioner(const bool& rebalance,
const std::string& reuse)
: prec_mp(Teuchos::null),
coords_mp(Teuchos::null),
rebalance_m(rebalance)
MueLuPreconditioner<Level>::MueLuPreconditioner(const std::string& reuse)
: prec_mp(Teuchos::null)
{
this->init_m(reuse);
}
......@@ -36,49 +33,11 @@ MueLuPreconditioner<Level>::MueLuPreconditioner(const bool& rebalance,
template <class Level>
void MueLuPreconditioner<Level>::create(const Teuchos::RCP<amr::matrix_t>& A,
Level* level_p)
Level* /*level_p*/)
{
coords_mp = Teuchos::null;
if ( rebalance_m ) {
const scalar_t* domain = level_p->geom.ProbLo();
const scalar_t* dx = level_p->cellSize();
coords_mp = Teuchos::rcp(
new amr::multivector_t(A->getDomainMap(), AMREX_SPACEDIM, false)
);
for (amrex::MFIter mfi(level_p->grids, level_p->dmap, true);
mfi.isValid(); ++mfi)
{
const AmrBox_t& tbx = mfi.tilebox();
const lo_t* lo = tbx.loVect();
const lo_t* hi = tbx.hiVect();
for (lo_t i = lo[0]; i <= hi[0]; ++i) {
for (lo_t j = lo[1]; j <= hi[1]; ++j) {
#if AMREX_SPACEDIM == 3
for (lo_t k = lo[2]; k <= hi[2]; ++k) {
#endif
AmrIntVect_t iv(D_DECL(i, j, k));
go_t gidx = level_p->serialize(iv);
coords_mp->replaceGlobalValue(gidx, 0, domain[0] + (0.5 + i) * dx[0]);
coords_mp->replaceGlobalValue(gidx, 1, domain[1] + (0.5 + j) * dx[1]);
#if AMREX_SPACEDIM == 3
coords_mp->replaceGlobalValue(gidx, 2, domain[2] + (0.5 + k) * dx[2]);
}
#endif
}
}
}
}
prec_mp = MueLu::CreateTpetraPreconditioner(A, params_m, coords_mp);
coords_mp = Teuchos::null;
typedef Tpetra::Operator<scalar_t, lo_t, go_t, amr::node_t> TpetraOperator_t;
Teuchos::RCP<TpetraOperator_t> At = Teuchos::rcp_dynamic_cast<TpetraOperator_t>(A);
prec_mp = MueLu::CreateTpetraPreconditioner(At, params_m);
}
......@@ -129,8 +88,8 @@ void MueLuPreconditioner<Level>::init_m(const std::string& reuse) {
params_m.set("sa: use filtered matrix", true);
params_m.set("filtered matrix: reuse eigenvalue", false); // false: more expensive
params_m.set("repartition: enable", rebalance_m);
params_m.set("repartition: rebalance P and R", rebalance_m);
params_m.set("repartition: enable", false);
params_m.set("repartition: rebalance P and R", false);
params_m.set("repartition: partitioner", "zoltan2");
params_m.set("repartition: min rows per proc", 800);
params_m.set("repartition: start level", 2);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment