Skip to content

GitLab

M
MoReFEM
Commit 60f3714a authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files
Options

#1535 Replace more Internal() by InternalForReadOnly().

parent 6e19b300
Hide whitespace changes
Inline Side-by-side
Showing with 40 additions and 29 deletions
Sources/ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hpp
View file @ 60f3714a
......@@ -156,14 +156,14 @@ namespace MoReFEM
*/
template
<
class MatrixT,
class MatrixU
class MatrixT,
class MatrixU
>
std::enable_if_t
<
std::is_base_of<Internal::Wrappers::Petsc::BaseMatrix, MatrixT>::value
&& std::is_base_of<Internal::Wrappers::Petsc::BaseMatrix, MatrixU>::value,
void
std::is_base_of<Internal::Wrappers::Petsc::BaseMatrix, MatrixT>::value
&& std::is_base_of<Internal::Wrappers::Petsc::BaseMatrix, MatrixU>::value,
void
>
MatTranspose(MatrixT& matrix1,
MatrixU& matrix2,
......@@ -403,6 +403,9 @@ namespace MoReFEM
* \param[in] matrix2 B in C = A * B^T.
* \param[out] matrix3 C in B in C = A * B^T. The matrix must be not allocated when this function is called.
* \copydoc doxygen_hide_invoking_file_and_line
*
* \attention It is not advised to use this wrapper (at least currently: as of PETSc v1.13.2, there is no support for this operation
* for MPIAIJ and MPIAIJ matrices - so the code you would write with it will not run in parallel.
*/
template
<
......
Sources/ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hxx
View file @ 60f3714a
......@@ -154,9 +154,9 @@ namespace MoReFEM
case DoReuseMatrix::yes:
{
result = out.Internal();
error_code = ::MatMatMatMult(matrix1.Internal(),
matrix2.Internal(),
matrix3.Internal(),
error_code = ::MatMatMatMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
matrix3.InternalForReadOnly(),
MAT_REUSE_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -166,9 +166,9 @@ namespace MoReFEM
}
case DoReuseMatrix::no:
{
error_code = ::MatMatMatMult(matrix1.Internal(),
matrix2.Internal(),
matrix3.Internal(),
error_code = ::MatMatMatMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
matrix3.InternalForReadOnly(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -261,29 +261,37 @@ namespace MoReFEM
const char* invoking_file, int invoking_line,
DoReuseMatrix do_reuse_matrix)
{
static_assert(std::is_const<MatrixU>() == false);
Mat result;
int error_code {};
switch (do_reuse_matrix)
{
case DoReuseMatrix::no:
{
error_code = ::MatTranspose(matrix1.Internal(), MAT_INITIAL_MATRIX, &result);
error_code = ::MatTranspose(matrix1.InternalForReadOnly(), MAT_INITIAL_MATRIX, &result);
matrix2.SetFromPetscMat(result);
break;
}
case DoReuseMatrix::yes:
{
result = matrix2.Internal();
error_code = ::MatTranspose(matrix1.Internal(), MAT_REUSE_MATRIX, &result);
error_code = ::MatTranspose(matrix1.InternalForReadOnly(), MAT_REUSE_MATRIX, &result);
break;
}
case DoReuseMatrix::in_place:
{
result = matrix1.Internal();
assert(matrix1.Internal() == matrix2.Internal() && "For in place transpose both arguments"
"are expected to be pointers to the"
"same PETSc matrix object.");
error_code = ::MatTranspose(matrix1.Internal(), MAT_INPLACE_MATRIX, &result);
result = matrix2.Internal();
assert(matrix1.InternalForReadOnly() == matrix2.Internal() && "For in place transpose both arguments"
"are expected to be pointers to the"
"same PETSc matrix object.");
// < note: Internal() and InternalForReadOnly() actually provide the same to the underlying
// pointer, but matrix1 might be const in some calls for the other values of do_reuse_matrix
// so Internal() couldn't be called on it reliably - to do so I would have to introduce
// a template class to enable partial specialization (I NEVER want to have to specify explicitly
// the matrix template parameters...).
error_code = ::MatTranspose(matrix2.Internal(), MAT_INPLACE_MATRIX, &result);
break;
}
} // switch
......@@ -368,8 +376,8 @@ namespace MoReFEM
{
result = matrix3.Internal();
error_code = ::MatTransposeMatMult(matrix1.Internal(),
matrix2.Internal(),
error_code = ::MatTransposeMatMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
MAT_REUSE_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -378,8 +386,8 @@ namespace MoReFEM
}
case DoReuseMatrix::no:
{
error_code = ::MatTransposeMatMult(matrix1.Internal(),
matrix2.Internal(),
error_code = ::MatTransposeMatMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -429,8 +437,8 @@ namespace MoReFEM
case DoReuseMatrix::yes:
{
result = matrix3.Internal();
error_code = ::MatMatTransposeMult(matrix1.Internal(),
matrix2.Internal(),
error_code = ::MatMatTransposeMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
MAT_REUSE_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -438,8 +446,8 @@ namespace MoReFEM
}
case DoReuseMatrix::no:
{
error_code = ::MatMatTransposeMult(matrix1.Internal(),
matrix2.Internal(),
error_code = ::MatMatTransposeMult(matrix1.InternalForReadOnly(),
matrix2.InternalForReadOnly(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
......@@ -528,7 +536,7 @@ namespace MoReFEM
IS *cperm,
const char* invoking_file, int invoking_line)
{
int error_code = ::MatGetOrdering(A.Internal(),
int error_code = ::MatGetOrdering(A.InternalForReadOnly(),
type,
rperm,
cperm);
......@@ -573,8 +581,8 @@ namespace MoReFEM
MatrixV& X,
const char* invoking_file, int invoking_line)
{
int error_code = ::MatMatSolve(A.Internal(),
B.Internal(),
int error_code = ::MatMatSolve(A.InternalForReadOnly(),
B.InternalForReadOnly(),
X.Internal());
if (error_code)
......
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