Commit 4e22d08b authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#723 ThirdParty/Matrix: define Matrix operations as free functions that act...

#723 ThirdParty/Matrix: define Matrix operations as free functions that act upon a template parameter of type Private::BaseMatrix. Hence ShellMatrix will be able to use them as well without creating any conflict with overloads that act upon vectors.
parent 041db26b
......@@ -339,6 +339,9 @@
BE55781C1BE250670097FB58 /* Solver.cpp in Sources */ = {isa = PBXBuildFile; fileRef = BE5578191BE250670097FB58 /* Solver.cpp */; };
BE55781D1BE250670097FB58 /* Solver.hpp in Headers */ = {isa = PBXBuildFile; fileRef = BE55781A1BE250670097FB58 /* Solver.hpp */; };
BE55781E1BE250670097FB58 /* Solver.hxx in Headers */ = {isa = PBXBuildFile; fileRef = BE55781B1BE250670097FB58 /* Solver.hxx */; };
BE56978E1BE383AD00B2EC67 /* MatrixOperations.hpp in Headers */ = {isa = PBXBuildFile; fileRef = BE56978B1BE383AD00B2EC67 /* MatrixOperations.hpp */; };
BE56978F1BE383AD00B2EC67 /* MatrixOperations.hxx in Headers */ = {isa = PBXBuildFile; fileRef = BE56978C1BE383AD00B2EC67 /* MatrixOperations.hxx */; };
BE5697A81BE3878700B2EC67 /* BaseMatrix.hpp in Headers */ = {isa = PBXBuildFile; fileRef = BE5697A51BE3878700B2EC67 /* BaseMatrix.hpp */; };
BE58EB751B9087B1006899EA /* InitializeHelper.cpp in Sources */ = {isa = PBXBuildFile; fileRef = BE58EB741B9087B1006899EA /* InitializeHelper.cpp */; };
BE58EB781B908827006899EA /* libModel.a in Frameworks */ = {isa = PBXBuildFile; fileRef = BE58EB4D1B908757006899EA /* libModel.a */; };
BE58EB7B1B90883E006899EA /* libModel.a in Frameworks */ = {isa = PBXBuildFile; fileRef = BE58EB4D1B908757006899EA /* libModel.a */; };
......@@ -3748,6 +3751,9 @@
BE5578191BE250670097FB58 /* Solver.cpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.cpp; name = Solver.cpp; path = Solver/Solver.cpp; sourceTree = "<group>"; };
BE55781A1BE250670097FB58 /* Solver.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = Solver.hpp; path = Solver/Solver.hpp; sourceTree = "<group>"; };
BE55781B1BE250670097FB58 /* Solver.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = Solver.hxx; path = Solver/Solver.hxx; sourceTree = "<group>"; };
BE56978B1BE383AD00B2EC67 /* MatrixOperations.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = MatrixOperations.hpp; sourceTree = "<group>"; };
BE56978C1BE383AD00B2EC67 /* MatrixOperations.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; path = MatrixOperations.hxx; sourceTree = "<group>"; };
BE5697A51BE3878700B2EC67 /* BaseMatrix.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = BaseMatrix.hpp; path = Private/BaseMatrix.hpp; sourceTree = "<group>"; };
BE589C281A160ABD00D23130 /* UniqueId.hpp */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = UniqueId.hpp; path = Sources/Utilities/UniqueId/UniqueId.hpp; sourceTree = "<group>"; };
BE589C291A160ABD00D23130 /* UniqueId.hxx */ = {isa = PBXFileReference; fileEncoding = 4; lastKnownFileType = sourcecode.cpp.h; name = UniqueId.hxx; path = Sources/Utilities/UniqueId/UniqueId.hxx; sourceTree = "<group>"; };
BE589C381A446A1200745D08 /* SConscript */ = {isa = PBXFileReference; explicitFileType = text.script.python; fileEncoding = 4; path = SConscript; sourceTree = "<group>"; };
......@@ -5760,6 +5766,14 @@
name = Solver;
sourceTree = "<group>";
};
BE56979D1BE3875D00B2EC67 /* Private */ = {
isa = PBXGroup;
children = (
BE5697A51BE3878700B2EC67 /* BaseMatrix.hpp */,
);
name = Private;
sourceTree = "<group>";
};
BE589C1E1A1609E100D23130 /* UniqueId */ = {
isa = PBXGroup;
children = (
......@@ -7992,6 +8006,9 @@
BEEFEF75196ECCC000C80FF1 /* MatrixPattern.cpp */,
BEEFEF76196ECCC000C80FF1 /* MatrixPattern.hpp */,
BEEFEF77196ECCC000C80FF1 /* MatrixPattern.hxx */,
BE56978B1BE383AD00B2EC67 /* MatrixOperations.hpp */,
BE56978C1BE383AD00B2EC67 /* MatrixOperations.hxx */,
BE56979D1BE3875D00B2EC67 /* Private */,
);
path = Matrix;
sourceTree = "<group>";
......@@ -8600,6 +8617,7 @@
BEC157F21A249C3D007D20EB /* Numeric.hpp in Headers */,
BE6E4EE31B2ABE8B0049BB2D /* AccessGhostContent.hxx in Headers */,
BE90E1A61A24929A00CCAFDE /* Helper.hxx in Headers */,
BE56978F1BE383AD00B2EC67 /* MatrixOperations.hxx in Headers */,
BE90E18F1A24929A00CCAFDE /* File.hpp in Headers */,
BE44C05B1AA463DF0030FA26 /* Pragma.hpp in Headers */,
BE4ED31D1A2CBAC400DE374E /* MatrixOperations.hpp in Headers */,
......@@ -8620,6 +8638,7 @@
BE90E1BB1A24929A00CCAFDE /* Snes.hpp in Headers */,
BE41A8CF1A24AA46004E4312 /* Mpi.hxx in Headers */,
BE1E87641B8DFB710002EE64 /* PrepareDefaultEntry.hxx in Headers */,
BE56978E1BE383AD00B2EC67 /* MatrixOperations.hpp in Headers */,
BE90E1861A24929A00CCAFDE /* ExceptionHelper.hxx in Headers */,
BE90E1CE1A2492AA00CCAFDE /* PetscSnes.hpp in Headers */,
BE41A8D21A24AA83004E4312 /* Petsc.hpp in Headers */,
......@@ -8667,6 +8686,7 @@
BE41A8CC1A24AA46004E4312 /* MpiScale.hpp in Headers */,
BE915A0C1AAF2E7000B4C474 /* LocalVectorStorage.hxx in Headers */,
BE90E1811A24929A00CCAFDE /* OpsFunction.hxx in Headers */,
BE5697A81BE3878700B2EC67 /* BaseMatrix.hpp in Headers */,
BE768E861B832795009B24CB /* Subtuple.hxx in Headers */,
BE768E7F1B8325BB009B24CB /* ManualParsing.hpp in Headers */,
BE6FFAD91A399E2200D048BD /* Numeric.hxx in Headers */,
......
......@@ -4,10 +4,12 @@
# include <vector>
# include "ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hpp"
# include "ThirdParty/Wrappers/Petsc/Snes.hpp"
# include "ThirdParty/Wrappers/Petsc/Print.hpp"
# include "ThirdParty/Wrappers/Petsc/Vector/AccessVectorContent.hpp"
# include "Utilities/Containers/UnorderedMap.hpp"
# include "Core/HappyHeart.hpp"
......
......@@ -183,6 +183,9 @@ namespace HappyHeart
void UpdateSolidDisplacementAfterSoF();
///@}
void dH();
private:
......
......@@ -46,28 +46,33 @@ namespace HappyHeart
UpdateSolidDisplacementAfterSoF();
auto& solid_displacement_after_dH = *solid_displacement_after_dH_;
solid_displacement_after_dH_->Copy(GetSolidDisplacementAfterSoF(), __FILE__, __LINE__);
dH();
}
template<class SolidVariationalFormulationPolicyT>
void Model<SolidVariationalFormulationPolicyT>::dH()
{
auto& solid_displacement_after_dH = *solid_displacement_after_dH_;
UpdateSolidVelocityOnInterfaceAfterDiffFluidImplicitStep(GetSolidDisplacementEvolutionInSoFOnInterface());
UpdateFluidDirichletCondition();
ImplicitFluidStep(differential::yes);
auto& fluid_varf = GetNonCstImplicitStepFluidVariationalFormulation();
{
const auto& fluid_sol = fluid_varf.GetSystemSolution(fluid_varf.GetNumberingSubset());
const PetscReal evaluation_state_min = fluid_sol.Min(__FILE__, __LINE__).second;
const PetscReal evaluation_state_max = fluid_sol.Max(__FILE__, __LINE__).second;
std::cout << "deltaVel EXTREMA ARE " << evaluation_state_min << " and " << evaluation_state_max
<< " on " << fluid_sol.GetProgramWiseSize(__FILE__, __LINE__) << " elements." << std::endl; // OK
}
const auto& differential_fluid_residual = fluid_varf.ComputeResidual();
......@@ -81,7 +86,7 @@ namespace HappyHeart
const auto& solid_varf = SolidVariationalFormulationPolicyT::GetVariationalFormulation();
const auto& solid_numbering_subset = solid_varf.GetNumberingSubset();
auto& solution = *solid_displacement_inside_dh_;
......@@ -101,12 +106,12 @@ namespace HappyHeart
const auto& solid_tangent = solid_varf.GetSystemMatrix(solid_numbering_subset, solid_numbering_subset);
// solid_tangent.View(mpi,__FILE__, __LINE__);
// solid_tangent.View(mpi,__FILE__, __LINE__);
const auto& god_of_dof = parent::GetGodOfDof(EnumUnderlyingType(MeshIndex::mesh));
auto& dirichlet_bc =
Private::DirichletBoundaryConditionManager::GetInstance().GetNonCstDirichletBoundaryCondition(EnumUnderlyingType(BoundaryConditionIndex::solid));
Private::DirichletBoundaryConditionManager::GetInstance().GetNonCstDirichletBoundaryCondition(EnumUnderlyingType(BoundaryConditionIndex::solid));
god_of_dof.ApplyBoundaryCondition(dirichlet_bc, solid_residual);
......@@ -118,7 +123,7 @@ namespace HappyHeart
{
const PetscReal evaluation_state_min = solution.Min(__FILE__, __LINE__).second;
const PetscReal evaluation_state_max = solution.Max(__FILE__, __LINE__).second;
std::cout << "diff solid solution EXTREMA ARE " << evaluation_state_min << " and " << evaluation_state_max
<< " on " << solution.GetProgramWiseSize(__FILE__, __LINE__) << " elements." << std::endl; // opposite of the Freefem output
}
......@@ -136,6 +141,7 @@ namespace HappyHeart
<< " on " << solid_displacement_after_dH.GetProgramWiseSize(__FILE__, __LINE__) << " elements." << std::endl; // opposite of the Freefem output
}
}
......
......@@ -6,6 +6,8 @@
// Copyright (c) 2014 Inria. All rights reserved.
//
#include "ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hpp"
#include "ModelInstances/Hyperelasticity/TimeSchemes/HalfSum/HalfSum.hpp"
......
......@@ -6,6 +6,8 @@
// Copyright (c) 2014 Inria. All rights reserved.
//
#include "ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hpp"
#include "ModelInstances/Hyperelasticity/TimeSchemes/Midpoint/Midpoint.hpp"
......
......@@ -315,34 +315,7 @@ namespace HappyHeart
Copy(source, invoking_file, invoking_line, SAME_NONZERO_PATTERN);
}
void AXPY(PetscScalar a, const Matrix& X, Matrix& Y, const char* invoking_file,
int invoking_line, const MatStructure& structure)
{
int error_code = ::MatAXPY(Y.Internal(), a, X.Internal(), structure);
if (error_code)
throw ExceptionNS::Exception(error_code, "MatAXPY", invoking_file, invoking_line);
}
void MatMult(const Matrix& matrix, const Vector& v1, Vector& v2, const char* invoking_file, int invoking_line)
{
int error_code = ::MatMult(matrix.Internal(), v1.Internal(), v2.Internal());
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMult", invoking_file, invoking_line);
}
void MatMultAdd(const Matrix& matrix, const Vector& v1,
const Vector& v2,
Vector& v3, const char* invoking_file, int invoking_line)
{
int error_code = ::MatMultAdd(matrix.Internal(), v1.Internal(), v2.Internal(), v3.Internal());
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMultAdd", invoking_file, invoking_line);
}
void Matrix::Scale(PetscScalar a, const char* invoking_file, int invoking_line)
{
......@@ -435,45 +408,6 @@ namespace HappyHeart
}
void MatMultTranspose(const Matrix& matrix, const Vector& v1, Vector& v2,
const char* invoking_file, int invoking_line)
{
int error_code = ::MatMultTranspose(matrix.Internal(), v1.Internal(), v2.Internal());
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMultTranspose", invoking_file, invoking_line);
}
void MatMultTransposeAdd(const Matrix& matrix, const Vector& v1, const Vector& v2, Vector& v3,
const char* invoking_file, int invoking_line)
{
int error_code = ::MatMultTransposeAdd(matrix.Internal(), v1.Internal(), v2.Internal(), v3.Internal());
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMultTransposeAdd", invoking_file, invoking_line);
}
void MatMatMult(const Matrix& matrix1,
const Matrix& matrix2,
Matrix& matrix3,
const char* invoking_file, int invoking_line)
{
Mat result;
int error_code = ::MatMatMult(matrix1.Internal(),
matrix2.Internal(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
matrix3.SetFromPetscMat(result);
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMatMult", invoking_file, invoking_line);
}
void Matrix::SetFromPetscMat(Mat petsc_matrix)
{
......@@ -482,32 +416,7 @@ namespace HappyHeart
petsc_matrix_ = petsc_matrix;
}
void MatMatMatMult(const Matrix& matrix1,
const Matrix& matrix2,
const Matrix& matrix3,
Matrix& matrix4,
const char* invoking_file, int invoking_line)
{
Mat result;
int error_code = ::MatMatMatMult(matrix1.Internal(),
matrix2.Internal(),
matrix3.Internal(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
matrix4.SetFromPetscMat(result);
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMatMatMult", invoking_file, invoking_line);
}
} // namespace Petsc
......
......@@ -18,6 +18,8 @@
# include "ThirdParty/IncludeWithoutWarning/Petsc/PetscSys.hpp"
# include "ThirdParty/IncludeWithoutWarning/Petsc/PetscMat.hpp"
# include "ThirdParty/Wrappers/Petsc/Matrix/Private/BaseMatrix.hpp"
# include "Utilities/SparseMatrix/CSRPattern.hpp"
......@@ -86,7 +88,7 @@ namespace HappyHeart
* by a 'Init***()' method or by using 'DuplicateLayout', 'Copy' or 'CompleteCopy'methods.
*
*/
class Matrix
class Matrix : private Private::BaseMatrix
{
public:
......@@ -150,6 +152,13 @@ namespace HappyHeart
void InitShellMatrix(unsigned int Nlocal_row, unsigned int Nlocal_column,
unsigned int Nglobal_row, unsigned int Nglobal_column
);
/*!
* \brief Handle over the internal Mat object.
*
......@@ -415,14 +424,13 @@ namespace HappyHeart
// \attention Do not forget to update Swap() if a new data member is added!
// =============================================================================
//! Underlying Petsc vector.
//! Underlying Petsc matrix.
Mat petsc_matrix_;
/*!
* \brief Whether the underlying Petsc matrix will be destroyed upon destruction of the object.
*
* Default behaviour is to do so, but in some cases (for instance when vector has been built from
* an existing Petsc Mat) it is wiser not to.
* Default behaviour is to do so, but in some cases it is wiser not to.
*/
bool do_petsc_destroy_;
......@@ -436,96 +444,7 @@ namespace HappyHeart
*/
void Swap(Matrix& A, Matrix& B);
/*!
* \brief Wrapper over MatAXPY, that performs Y = a * X + Y.
*
* \param[in] structure either SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN or SUBSET_NONZERO_PATTERN.
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void AXPY(PetscScalar a, const Matrix& X, Matrix& Y, const char* invoking_file, int invoking_line,
const MatStructure& structure = SAME_NONZERO_PATTERN);
/*!
* \brief Wrapper over MatMult, that performs v2 = matrix * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMult(const Matrix& matrix, const Vector& v1, Vector& v2,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultAdd, that performs v3 = v2 + matrix * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMultAdd(const Matrix& matrix, const Vector& v1, const Vector& v2, Vector& v3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultTranspose, that performs v2 = transpose(matrix) * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMultTranspose(const Matrix& matrix, const Vector& v1, Vector& v2,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultTransposeAdd, that performs v3 = v2 + transpose(matrix) * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMultTransposeAdd(const Matrix& matrix, const Vector& v1, const Vector& v2, Vector& v3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMatMult, that performs C = A * B.
*
* \todo #684 Investigate to use the argument fill, which provides an estimation of the non zero of the
* resulting matrix. Currently PETSC_DEFAULT is used.
* One can also reuse a pattern for several matrices, but so far I do not need that in HappyHeart so
* I equally use the default setting.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMatMult(const Matrix& matrix1,
const Matrix& matrix2,
Matrix& matrix3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMatMatMult, that performs D = A * B * C.
*
* \todo #684 Investigate to use the argument fill, which provides an estimation of the non zero of the
* resulting matrix. Currently PETSC_DEFAULT is used.
* One can also reuse a pattern for several matrices, but so far I do not need that in HappyHeart so
* I equally use the default setting.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
void MatMatMatMult(const Matrix& matrix1,
const Matrix& matrix2,
const Matrix& matrix3,
Matrix& matrix4,
const char* invoking_file, int invoking_line);
} // namespace Petsc
......
......@@ -31,6 +31,8 @@ namespace HappyHeart
}
} // namespace Petsc
......
//
// MatrixTOperations.hpp
// HappyHeart
//
// Created by Sebastien Gilles on 30/10/15.
// Copyright © 2015 Inria. All rights reserved.
//
#ifndef HAPPY_HEART_x_THIRD_PARTY_x_WRAPPERS_x_PETSC_x_MATRIX_x_MATRIX_OPERATIONS_HPP_
# define HAPPY_HEART_x_THIRD_PARTY_x_WRAPPERS_x_PETSC_x_MATRIX_x_MATRIX_OPERATIONS_HPP_
# include <memory>
# include <vector>
# include "ThirdParty/IncludeWithoutWarning/Petsc/PetscSys.hpp"
# include "ThirdParty/IncludeWithoutWarning/Petsc/PetscMat.hpp"
# include "ThirdParty/Wrappers/Petsc/Vector/Vector.hpp"
# include "ThirdParty/Wrappers/Petsc/Matrix/Private/BaseMatrix.hpp"
# include "ThirdParty/Wrappers/Petsc/Exceptions/Petsc.hpp"
namespace HappyHeart
{
namespace Wrappers
{
namespace Petsc
{
/*!
* \brief Wrapper over MatAXPY, that performs Y = a * X + Y.
*
* \param[in] structure either SAME_NONZERO_PATTERN, DIFFERENT_NONZERO_PATTERN or SUBSET_NONZERO_PATTERN.
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
AXPY(PetscScalar a,
const MatrixT& X,
MatrixT& Y,
const char* invoking_file, int invoking_line,
const MatStructure& structure = SAME_NONZERO_PATTERN);
/*!
* \brief Wrapper over MatMult, that performs v2 = matrix * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMult(const MatrixT& matrix, const Vector& v1, Vector& v2,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultAdd, that performs v3 = v2 + matrix * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMultAdd(const MatrixT& matrix, const Vector& v1, const Vector& v2, Vector& v3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultTranspose, that performs v2 = transpose(matrix) * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMultTranspose(const MatrixT& matrix, const Vector& v1, Vector& v2,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMultTransposeAdd, that performs v3 = v2 + transpose(matrix) * v1.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMultTransposeAdd(const MatrixT& matrix, const Vector& v1, const Vector& v2, Vector& v3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMatMult, that performs C = A * B.
*
* \todo #684 Investigate to use the argument fill, which provides an estimation of the non zero of the
* resulting matrix. Currently PETSC_DEFAULT is used.
* One can also reuse a pattern for several matrices, but so far I do not need that in HappyHeart so
* I equally use the default setting.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMatMult(const MatrixT& matrix1,
const MatrixT& matrix2,
MatrixT& matrix3,
const char* invoking_file, int invoking_line);
/*!
* \brief Wrapper over MatMatMatMult, that performs D = A * B * C.
*
* \todo #684 Investigate to use the argument fill, which provides an estimation of the non zero of the
* resulting matrix. Currently PETSC_DEFAULT is used.
* One can also reuse a pattern for several matrices, but so far I do not need that in HappyHeart so
* I equally use the default setting.
*
* \param[in] invoking_file File that invoked the function or class; usually __FILE__.
* \param[in] invoking_line File that invoked the function or class; usually __LINE__.
*/
template<class MatrixT>
std::enable_if_t<std::is_base_of<Private::BaseMatrix, MatrixT>::value, void>
MatMatMatMult(const MatrixT& matrix1,
const MatrixT& matrix2,
const MatrixT& matrix3,
MatrixT& matrix4,
const char* invoking_file, int invoking_line);