Skip to content
Merged
Show file tree
Hide file tree
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
Original file line number Diff line number Diff line change
Expand Up @@ -7,15 +7,15 @@
#include "dmumps_c.h"
#include "mpi.h"

class DirectSolverGive : public DirectSolver
class DirectSolver_COO_MUMPS_Give : public DirectSolver
{
public:
explicit DirectSolverGive(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients, bool DirBC_Interior,
int num_omp_threads);
explicit DirectSolver_COO_MUMPS_Give(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads);

~DirectSolverGive() override;
~DirectSolver_COO_MUMPS_Give() override;
// Note: The rhs (right-hand side) vector gets overwritten during the solution process.
void solveInPlace(Vector<double> solution) override;

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -7,15 +7,15 @@
#include "dmumps_c.h"
#include "mpi.h"

class DirectSolverTake : public DirectSolver
class DirectSolver_COO_MUMPS_Take : public DirectSolver
{
public:
explicit DirectSolverTake(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients, bool DirBC_Interior,
int num_omp_threads);
explicit DirectSolver_COO_MUMPS_Take(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads);

~DirectSolverTake() override;
~DirectSolver_COO_MUMPS_Take() override;
// Note: The rhs (right-hand side) vector gets overwritten during the solution process.
void solveInPlace(Vector<double> solution) override;

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -3,15 +3,15 @@
#include "../directSolver.h"
#include <Kokkos_Core.hpp>

class DirectSolverGiveCustomLU : public DirectSolver
class DirectSolver_CSR_LU_Give : public DirectSolver
{
public:
explicit DirectSolverGiveCustomLU(const PolarGrid& grid, const LevelCache& level_cache,
explicit DirectSolver_CSR_LU_Give(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads);

~DirectSolverGiveCustomLU() override;
~DirectSolver_CSR_LU_Give() override;
// Note: The rhs (right-hand side) vector gets overwritten with the solution.
void solveInPlace(Vector<double> solution) override;

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,15 +2,15 @@

#include "../directSolver.h"

class DirectSolverTakeCustomLU : public DirectSolver
class DirectSolver_CSR_LU_Take : public DirectSolver
{
public:
explicit DirectSolverTakeCustomLU(const PolarGrid& grid, const LevelCache& level_cache,
explicit DirectSolver_CSR_LU_Take(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads);

~DirectSolverTakeCustomLU() override;
~DirectSolver_CSR_LU_Take() override;
// Note: The rhs (right-hand side) vector gets overwritten with the solution.
void solveInPlace(Vector<double> solution) override;

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
/* Boundary Symmetry Shift */
/* ----------------------- */

void DirectSolverGive::applySymmetryShiftInnerBoundary(Vector<double> x) const
void DirectSolver_COO_MUMPS_Give::applySymmetryShiftInnerBoundary(Vector<double> x) const
{
assert(DirBC_Interior_);

Expand Down Expand Up @@ -67,7 +67,7 @@ void DirectSolverGive::applySymmetryShiftInnerBoundary(Vector<double> x) const
}
}

void DirectSolverGive::applySymmetryShiftOuterBoundary(Vector<double> x) const
void DirectSolver_COO_MUMPS_Give::applySymmetryShiftOuterBoundary(Vector<double> x) const
{
int i_r;
double r;
Expand Down Expand Up @@ -123,7 +123,7 @@ void DirectSolverGive::applySymmetryShiftOuterBoundary(Vector<double> x) const
}

// clang-format off
void DirectSolverGive::applySymmetryShift(Vector<double> x) const
void DirectSolver_COO_MUMPS_Give::applySymmetryShift(Vector<double> x) const
{
assert(std::ssize(x) == grid_.numberOfNodes());
assert(grid_.nr() >= 4);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -765,7 +765,7 @@
} \
} while (0)

void DirectSolverGive::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Give::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCOO<double>& solver_matrix)
{
const double r = grid_.radius(i_r);
for (int i_theta = 0; i_theta < grid_.ntheta(); i_theta++) {
Expand All @@ -781,7 +781,8 @@ void DirectSolverGive::buildSolverMatrixCircleSection(const int i_r, SparseMatri
}
}

void DirectSolverGive::buildSolverMatrixRadialSection(const int i_theta, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Give::buildSolverMatrixRadialSection(const int i_theta,
SparseMatrixCOO<double>& solver_matrix)
{
const double theta = grid_.theta(i_theta);
for (int i_r = grid_.numberSmootherCircles(); i_r < grid_.nr(); i_r++) {
Expand All @@ -801,7 +802,7 @@ void DirectSolverGive::buildSolverMatrixRadialSection(const int i_theta, SparseM

/* ------------------------------------------------------------------------ */
/* If the indexing is not smoother-based, please adjust the access patterns */
SparseMatrixCOO<double> DirectSolverGive::buildSolverMatrix()
SparseMatrixCOO<double> DirectSolver_COO_MUMPS_Give::buildSolverMatrix()
{
const int n = grid_.numberOfNodes();
const int nnz = getNonZeroCountSolverMatrix();
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,17 +2,17 @@

#ifdef GMGPOLAR_USE_MUMPS

DirectSolverGive::DirectSolverGive(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients, bool DirBC_Interior,
int num_omp_threads)
DirectSolver_COO_MUMPS_Give::DirectSolver_COO_MUMPS_Give(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads)
: DirectSolver(grid, level_cache, domain_geometry, density_profile_coefficients, DirBC_Interior, num_omp_threads)
{
solver_matrix_ = buildSolverMatrix();
initializeMumpsSolver(mumps_solver_, solver_matrix_);
}

void DirectSolverGive::solveInPlace(Vector<double> solution)
void DirectSolver_COO_MUMPS_Give::solveInPlace(Vector<double> solution)
{
// Adjusts the right-hand side vector to account for symmetry corrections.
// This transforms the system matrixA * solution = rhs into the equivalent system:
Expand All @@ -24,7 +24,7 @@ void DirectSolverGive::solveInPlace(Vector<double> solution)
solveWithMumps(solution);
}

DirectSolverGive::~DirectSolverGive()
DirectSolver_COO_MUMPS_Give::~DirectSolver_COO_MUMPS_Give()
{
finalizeMumpsSolver(mumps_solver_);
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,8 @@

#ifdef GMGPOLAR_USE_MUMPS

void DirectSolverGive::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Give::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver,
SparseMatrixCOO<double>& solver_matrix)
{
/*
* MUMPS (a parallel direct solver) uses 1-based indexing,
Expand Down Expand Up @@ -102,7 +103,7 @@ void DirectSolverGive::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver, Spars
}
}

void DirectSolverGive::solveWithMumps(Vector<double> result_rhs)
void DirectSolver_COO_MUMPS_Give::solveWithMumps(Vector<double> result_rhs)
{
mumps_solver_.job = JOB_COMPUTE_SOLUTION;
mumps_solver_.nrhs = 1;
Expand All @@ -115,7 +116,7 @@ void DirectSolverGive::solveWithMumps(Vector<double> result_rhs)
}
}

void DirectSolverGive::finalizeMumpsSolver(DMUMPS_STRUC_C& mumps_solver)
void DirectSolver_COO_MUMPS_Give::finalizeMumpsSolver(DMUMPS_STRUC_C& mumps_solver)
{
mumps_solver.job = JOB_END;
dmumps_c(&mumps_solver);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

#ifdef GMGPOLAR_USE_MUMPS

const Stencil& DirectSolverGive::getStencil(int i_r) const
const Stencil& DirectSolver_COO_MUMPS_Give::getStencil(int i_r) const
{
assert(0 <= i_r && i_r < grid_.nr());
assert(grid_.nr() >= 4);
Expand All @@ -25,7 +25,7 @@ const Stencil& DirectSolverGive::getStencil(int i_r) const
throw std::out_of_range("Invalid index for stencil");
}

int DirectSolverGive::getNonZeroCountSolverMatrix() const
int DirectSolver_COO_MUMPS_Give::getNonZeroCountSolverMatrix() const
{
const int size_stencil_inner_boundary = DirBC_Interior_ ? 1 : 7;
const int size_stencil_next_inner_boundary = DirBC_Interior_ ? 6 : 9;
Expand All @@ -42,7 +42,7 @@ int DirectSolverGive::getNonZeroCountSolverMatrix() const

/* ----------------------------------------------------------------- */
/* If the indexing is not smoother-based, please adjust the indexing */
int DirectSolverGive::getSolverMatrixIndex(const int i_r, const int i_theta) const
int DirectSolver_COO_MUMPS_Give::getSolverMatrixIndex(const int i_r, const int i_theta) const
{
const int size_stencil_inner_boundary = DirBC_Interior_ ? 1 : 7;
const int size_stencil_next_inner_boundary = DirBC_Interior_ ? 6 : 9;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
/* Boundary Symmetry Shift */
/* ----------------------- */

void DirectSolverTake::applySymmetryShiftInnerBoundary(Vector<double> x) const
void DirectSolver_COO_MUMPS_Take::applySymmetryShiftInnerBoundary(Vector<double> x) const
{
assert(DirBC_Interior_);

Expand Down Expand Up @@ -45,7 +45,7 @@ void DirectSolverTake::applySymmetryShiftInnerBoundary(Vector<double> x) const
}
}

void DirectSolverTake::applySymmetryShiftOuterBoundary(Vector<double> x) const
void DirectSolver_COO_MUMPS_Take::applySymmetryShiftOuterBoundary(Vector<double> x) const
{
assert(level_cache_.cacheDensityProfileCoefficients());
assert(level_cache_.cacheDomainGeometry());
Expand Down Expand Up @@ -80,7 +80,7 @@ void DirectSolverTake::applySymmetryShiftOuterBoundary(Vector<double> x) const
}
}

void DirectSolverTake::applySymmetryShift(Vector<double> x) const
void DirectSolver_COO_MUMPS_Take::applySymmetryShift(Vector<double> x) const
{
assert(std::ssize(x) == grid_.numberOfNodes());
assert(grid_.nr() >= 4);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -447,7 +447,7 @@
} \
} while (0)

void DirectSolverTake::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Take::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCOO<double>& solver_matrix)
{
assert(level_cache_.cacheDensityProfileCoefficients());
assert(level_cache_.cacheDomainGeometry());
Expand All @@ -465,7 +465,8 @@ void DirectSolverTake::buildSolverMatrixCircleSection(const int i_r, SparseMatri
}
}

void DirectSolverTake::buildSolverMatrixRadialSection(const int i_theta, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Take::buildSolverMatrixRadialSection(const int i_theta,
SparseMatrixCOO<double>& solver_matrix)
{
assert(level_cache_.cacheDensityProfileCoefficients());
assert(level_cache_.cacheDomainGeometry());
Expand All @@ -487,7 +488,7 @@ void DirectSolverTake::buildSolverMatrixRadialSection(const int i_theta, SparseM

/* ------------------------------------------------------------------------ */
/* If the indexing is not smoother-based, please adjust the access patterns */
SparseMatrixCOO<double> DirectSolverTake::buildSolverMatrix()
SparseMatrixCOO<double> DirectSolver_COO_MUMPS_Take::buildSolverMatrix()
{
const int n = grid_.numberOfNodes();
const int nnz = getNonZeroCountSolverMatrix();
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,17 +2,17 @@

#ifdef GMGPOLAR_USE_MUMPS

DirectSolverTake::DirectSolverTake(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients, bool DirBC_Interior,
int num_omp_threads)
DirectSolver_COO_MUMPS_Take::DirectSolver_COO_MUMPS_Take(const PolarGrid& grid, const LevelCache& level_cache,
const DomainGeometry& domain_geometry,
const DensityProfileCoefficients& density_profile_coefficients,
bool DirBC_Interior, int num_omp_threads)
: DirectSolver(grid, level_cache, domain_geometry, density_profile_coefficients, DirBC_Interior, num_omp_threads)
{
solver_matrix_ = buildSolverMatrix();
initializeMumpsSolver(mumps_solver_, solver_matrix_);
}

void DirectSolverTake::solveInPlace(Vector<double> solution)
void DirectSolver_COO_MUMPS_Take::solveInPlace(Vector<double> solution)
{
// Adjusts the right-hand side vector to account for symmetry corrections.
// This transforms the system matrixA * solution = rhs into the equivalent system:
Expand All @@ -24,7 +24,7 @@ void DirectSolverTake::solveInPlace(Vector<double> solution)
solveWithMumps(solution);
}

DirectSolverTake::~DirectSolverTake()
DirectSolver_COO_MUMPS_Take::~DirectSolver_COO_MUMPS_Take()
{
finalizeMumpsSolver(mumps_solver_);
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,8 @@

#ifdef GMGPOLAR_USE_MUMPS

void DirectSolverTake::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver, SparseMatrixCOO<double>& solver_matrix)
void DirectSolver_COO_MUMPS_Take::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver,
SparseMatrixCOO<double>& solver_matrix)
{
/*
* MUMPS (a parallel direct solver) uses 1-based indexing,
Expand Down Expand Up @@ -102,7 +103,7 @@ void DirectSolverTake::initializeMumpsSolver(DMUMPS_STRUC_C& mumps_solver, Spars
}
}

void DirectSolverTake::solveWithMumps(Vector<double> result_rhs)
void DirectSolver_COO_MUMPS_Take::solveWithMumps(Vector<double> result_rhs)
{
mumps_solver_.job = JOB_COMPUTE_SOLUTION;
mumps_solver_.nrhs = 1;
Expand All @@ -115,7 +116,7 @@ void DirectSolverTake::solveWithMumps(Vector<double> result_rhs)
}
}

void DirectSolverTake::finalizeMumpsSolver(DMUMPS_STRUC_C& mumps_solver)
void DirectSolver_COO_MUMPS_Take::finalizeMumpsSolver(DMUMPS_STRUC_C& mumps_solver)
{
mumps_solver.job = JOB_END;
dmumps_c(&mumps_solver);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

#ifdef GMGPOLAR_USE_MUMPS

const Stencil& DirectSolverTake::getStencil(int i_r) const
const Stencil& DirectSolver_COO_MUMPS_Take::getStencil(int i_r) const
{
assert(0 <= i_r && i_r < grid_.nr());
assert(grid_.nr() >= 4);
Expand All @@ -25,7 +25,7 @@ const Stencil& DirectSolverTake::getStencil(int i_r) const
throw std::out_of_range("Invalid index for stencil");
}

int DirectSolverTake::getNonZeroCountSolverMatrix() const
int DirectSolver_COO_MUMPS_Take::getNonZeroCountSolverMatrix() const
{
const int size_stencil_inner_boundary = DirBC_Interior_ ? 1 : 7;
const int size_stencil_next_inner_boundary = DirBC_Interior_ ? 6 : 9;
Expand All @@ -42,7 +42,7 @@ int DirectSolverTake::getNonZeroCountSolverMatrix() const

/* ----------------------------------------------------------------- */
/* If the indexing is not smoother-based, please adjust the indexing */
int DirectSolverTake::getSolverMatrixIndex(const int i_r, const int i_theta) const
int DirectSolver_COO_MUMPS_Take::getSolverMatrixIndex(const int i_r, const int i_theta) const
{
const int size_stencil_inner_boundary = DirBC_Interior_ ? 1 : 7;
const int size_stencil_next_inner_boundary = DirBC_Interior_ ? 6 : 9;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -730,7 +730,7 @@
} \
} while (0)

void DirectSolverGiveCustomLU::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCSR<double>& solver_matrix)
void DirectSolver_CSR_LU_Give::buildSolverMatrixCircleSection(const int i_r, SparseMatrixCSR<double>& solver_matrix)
{
const double r = grid_.radius(i_r);
for (int i_theta = 0; i_theta < grid_.ntheta(); i_theta++) {
Expand All @@ -746,7 +746,7 @@ void DirectSolverGiveCustomLU::buildSolverMatrixCircleSection(const int i_r, Spa
}
}

void DirectSolverGiveCustomLU::buildSolverMatrixRadialSection(const int i_theta, SparseMatrixCSR<double>& solver_matrix)
void DirectSolver_CSR_LU_Give::buildSolverMatrixRadialSection(const int i_theta, SparseMatrixCSR<double>& solver_matrix)
{
const double theta = grid_.theta(i_theta);
for (int i_r = grid_.numberSmootherCircles(); i_r < grid_.nr(); i_r++) {
Expand All @@ -766,7 +766,7 @@ void DirectSolverGiveCustomLU::buildSolverMatrixRadialSection(const int i_theta,

/* ------------------------------------------------------------------------ */
/* If the indexing is not smoother-based, please adjust the access patterns */
SparseMatrixCSR<double> DirectSolverGiveCustomLU::buildSolverMatrix()
SparseMatrixCSR<double> DirectSolver_CSR_LU_Give::buildSolverMatrix()
{
const int n = grid_.numberOfNodes();

Expand Down
Loading