Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
104 changes: 65 additions & 39 deletions examples/monodomain_DG2D.cc
Original file line number Diff line number Diff line change
Expand Up @@ -61,7 +61,7 @@ static constexpr unsigned int starting_level = 1;

// matrix-free related parameters
static constexpr bool use_matrix_free_action = true;
static constexpr unsigned int degree_finite_element = 4;
static constexpr unsigned int degree_finite_element = 1;
constexpr unsigned int n_qpoints = degree_finite_element + 1;
static constexpr unsigned int n_components = 1;

Expand All @@ -82,22 +82,8 @@ namespace Utils
initialize(const MatrixType &matrix)
{
coarse_matrix = &matrix;
if constexpr (std::is_same_v<PreconditionerType, SparseDirectMUMPS>)
{
#ifdef DEAL_II_WITH_MUMPS
precondition = std::make_unique<SparseDirectMUMPS>(
typename SparseDirectMUMPS::AdditionalData(),
matrix.get_mpi_communicator());
precondition->initialize(*coarse_matrix);
#else
AssertThrow(false,
ExcMessage(
"MUMPS support is required but not available. "
"Please recompile deal.II with MUMPS enabled."));
#endif
}
else if constexpr (std::is_same_v<PreconditionerType,
TrilinosWrappers::PreconditionAMG>)
if constexpr (std::is_same_v<PreconditionerType,
TrilinosWrappers::PreconditionAMG>)
{
precondition = std::make_unique<TrilinosWrappers::PreconditionAMG>();
precondition->initialize(*coarse_matrix);
Expand Down Expand Up @@ -126,23 +112,6 @@ namespace Utils
start = MPI_Wtime();
#endif
solver_coarse.solve(*coarse_matrix, dst, src, *precondition);
#ifdef AGGLO_DEBUG
stop = MPI_Wtime();
#endif
}
else if constexpr (std::is_same_v<PreconditionerType, SparseDirectMUMPS>)
{
#ifdef AGGLO_DEBUG
start = MPI_Wtime();
#endif
#ifdef DEAL_II_WITH_MUMPS
precondition->vmult(dst, src);
#else
AssertThrow(false,
ExcMessage(
"MUMPS support is required but not available. "
"Please recompile deal.II with MUMPS enabled."));
#endif
#ifdef AGGLO_DEBUG
stop = MPI_Wtime();
#endif
Expand All @@ -163,6 +132,63 @@ namespace Utils
std::unique_ptr<PreconditionerType> precondition;
};



template <typename MatrixType, typename Number, typename DirectSolverType>
class MGCoarseDirectMUMPS
: public MGCoarseGridBase<LinearAlgebra::distributed::Vector<Number>>
{
public:
MGCoarseDirectMUMPS()
{}

void
initialize(const MatrixType &matrix)
{
#ifdef DEAL_II_WITH_MUMPS
coarse_matrix = &matrix;
direct_solver = std::make_unique<SparseDirectMUMPS>(
typename SparseDirectMUMPS::AdditionalData(),
matrix.get_mpi_communicator());
direct_solver->initialize(*coarse_matrix);

#else
DEAL_II_NOT_IMPLEMENTED();
#endif
}

virtual void
operator()(
const unsigned int level,
LinearAlgebra::distributed::Vector<double> &dst,
const LinearAlgebra::distributed::Vector<double> &src) const override
{
#ifdef DEAL_II_WITH_MUMPS
# ifdef AGGLO_DEBUG
(void)level;
[[maybe_unused]] double start, stop;
start = MPI_Wtime();
# endif
direct_solver->vmult(dst, src);
# ifdef AGGLO_DEBUG
stop = MPI_Wtime();
# endif


# ifdef AGGLO_DEBUG
if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
std::cout << "Coarse solver elapsed time: " << stop - start << "[s]"
<< std::endl;
# endif
#else
DEAL_II_NOT_IMPLEMENTED();
#endif
}

const MatrixType *coarse_matrix;
std::unique_ptr<DirectSolverType> direct_solver;
};

} // namespace Utils

// Model parameters for FitzHugh-Nagumo
Expand Down Expand Up @@ -499,9 +525,9 @@ class IonicModel


// Multigrid related types
using LevelMatrixType = LinearOperatorMG<VectorType, VectorType>;
using SmootherType = PreconditionChebyshev<LevelMatrixType, VectorType>;
using CoarseSolverType = TrilinosWrappers::PreconditionAMG;
using LevelMatrixType = LinearOperatorMG<VectorType, VectorType>;
using SmootherType = PreconditionChebyshev<LevelMatrixType, VectorType>;
using CoarseSolver = TrilinosWrappers::PreconditionAMG;

MGLevelObject<std::unique_ptr<TrilinosWrappers::SparseMatrix>>
multigrid_matrices;
Expand All @@ -520,7 +546,7 @@ class IonicModel

std::unique_ptr<Utils::MGCoarseIterative<TrilinosWrappers::SparseMatrix,
double,
CoarseSolverType>>
CoarseSolver>>
mg_coarse;

std::unique_ptr<mg::SmootherRelaxation<SmootherType, VectorType>> mg_smoother;
Expand Down Expand Up @@ -894,7 +920,7 @@ IonicModel<dim>::setup_multigrid()
mg_coarse =
std::make_unique<Utils::MGCoarseIterative<TrilinosWrappers::SparseMatrix,
double,
CoarseSolverType>>();
CoarseSolver>>();
mg_coarse->initialize(*multigrid_matrices[min_level]);

pcout << "Coarse solver initialized" << std::endl;
Expand Down
Loading