MatrixVectorTraits.cpp

Go to the documentation of this file.

// SPDX-FileCopyrightText: Copyright (c) OpenGeoSys Community (opengeosys.org)
// SPDX-License-Identifier: BSD-3-Clause
 
#include "MatrixVectorTraits.h"
 
#include <cassert>
 
#include "MatrixSpecifications.h"
 
#ifdef USE_PETSC
 
namespace MathLib
{
std::unique_ptr<PETScMatrix> MatrixVectorTraits<PETScMatrix>::newInstance()
{
    return std::make_unique<PETScMatrix>();
}
 
std::unique_ptr<PETScMatrix> MatrixVectorTraits<PETScMatrix>::newInstance(
    PETScMatrix const& A)
{
    return std::make_unique<PETScMatrix>(A);
}
 
std::unique_ptr<PETScMatrix> MatrixVectorTraits<PETScMatrix>::newInstance(
    MatrixSpecifications const& spec)
{
    auto const nrows = spec.nrows;
    auto const ncols = spec.ncols;
 
    if (spec.sparsity_pattern)
    {
        // Assert that the misuse of the sparsity pattern is consistent.
        assert(spec.sparsity_pattern->size() == 1);
 
        auto const max_nonzeroes = spec.sparsity_pattern->front();
 
        PETScMatrixOption mat_opt;
        mat_opt.d_nz = max_nonzeroes;
        mat_opt.o_nz = max_nonzeroes;
        mat_opt.is_global_size = false;
        return std::make_unique<PETScMatrix>(nrows, ncols, mat_opt);
    }
    else
        return std::make_unique<PETScMatrix>(nrows, ncols);
}
 
std::unique_ptr<PETScVector> MatrixVectorTraits<PETScVector>::newInstance()
{
    return std::make_unique<PETScVector>();
}
 
std::unique_ptr<PETScVector> MatrixVectorTraits<PETScVector>::newInstance(
    PETScVector const& x)
{
    return std::make_unique<PETScVector>(x);
}
 
std::unique_ptr<PETScVector> MatrixVectorTraits<PETScVector>::newInstance(
    MatrixSpecifications const& spec)
{
    auto const is_global_size = false;
 
    if (spec.ghost_indices != nullptr)
    {
        return std::make_unique<PETScVector>(spec.nrows, *spec.ghost_indices,
                                             is_global_size);
    }
    else
    {
        return std::make_unique<PETScVector>(spec.nrows, is_global_size);
    }
}
 
std::unique_ptr<PETScVector> MatrixVectorTraits<PETScVector>::newInstance(
    PETScVector::IndexType const length)
{
    auto const is_global_size = false;
 
    return std::make_unique<PETScVector>(length, is_global_size);
}
}  // namespace MathLib
 
#else
 
namespace MathLib
{
std::unique_ptr<EigenMatrix> MatrixVectorTraits<EigenMatrix>::newInstance()
{
    return std::make_unique<EigenMatrix>(0, 0);  // TODO default constructor
}
 
std::unique_ptr<EigenMatrix> MatrixVectorTraits<EigenMatrix>::newInstance(
    EigenMatrix const& A)
{
    return std::make_unique<EigenMatrix>(A);
}
 
std::unique_ptr<EigenMatrix> MatrixVectorTraits<EigenMatrix>::newInstance(
    MatrixSpecifications const& spec)
{
    auto A = std::make_unique<EigenMatrix>(spec.nrows);
 
    if (spec.sparsity_pattern)
    {
        setMatrixSparsity(*A, *spec.sparsity_pattern);
    }
 
    return A;
}
 
std::unique_ptr<EigenVector> MatrixVectorTraits<EigenVector>::newInstance()
{
    return std::make_unique<EigenVector>();
}
 
std::unique_ptr<EigenVector> MatrixVectorTraits<EigenVector>::newInstance(
    EigenVector const& x)
{
    return std::make_unique<EigenVector>(x);
}
 
std::unique_ptr<EigenVector> MatrixVectorTraits<EigenVector>::newInstance(
    MatrixSpecifications const& spec)
{
    return std::make_unique<EigenVector>(spec.nrows);
}
 
std::unique_ptr<EigenVector> MatrixVectorTraits<EigenVector>::newInstance(
    Eigen::SparseMatrix<double>::Index const length)
{
    return std::make_unique<EigenVector>(length);
}
}  // namespace MathLib
 
#endif